在R中组合3个不同的矩阵图

Question

谁能告诉我如何创建一个具有3个不同矩阵数据集的图。一般来说，我有3个不同的数据所有1 * 1001尺寸的矩阵，我想在同一个图上绘制所有3个数据。

我设法得到一个矩阵来一次绘制，并组装代码来创建其他2个矩阵，但不能绘制它。 B [i，]是随机生成的数据。我想知道的是在一张图上将所有3个图组合在一起的代码是什么。

代码用于一个矩阵： n时间< -1000 average.price.at.each.timestep < - 矩阵（0，nrow = 1，NcoI位= n时间+ 1）

for(i in 1:(ntime+1)){ 
average.price.at.each.timestep[i]<-mean(B[i,]) 
} 

matplot(t, t(average.price.at.each.timestep), type="l", lty=1, main="MC Price of a Zero Coupon Bond", ylab="Price", xlab = "Option Exercise Date") 

代码3 ：

average.price.at.each.timestep<-matrix(0,nrow=1,ncol=ntime+1) 
s.e.at.each.time <-matrix(0,nrow=1,ncol=ntime+1) 
upper.c.l.at <- matrix(0,nrow=1,ncol=ntime+1) 
lower.c.l.at <- matrix(0,nrow=1,ncol=ntime+1) 
std <- function(x) sd(x)/sqrt(length(x)) 

for(i in 1:(ntime+1)){ 
average.price.at.each.timestep[i]<-mean(B[i,]) 
s.e.at.each.time[i] <- std(B[i,]) 
upper.c.l.at[i] <- average.price.at.each.timestep[i]+1.96*s.e.at.each.time[i] 
lower.c.l.at[i] <- average.price.at.each.timestep[i]-1.96*s.e.at.each.time[i] 
} 

我仍然为此而努力，因为我不能让给予配合我的数据集的解决方案，我现在已经包括下面的代码，作为一个工作示例生成矩阵B，所以你可以看到我正在处理的数据。正如你可以看到它产生了不同价格的情节，我想要一个平均价格和置信区间的平均值。

# Define Bond Price Parameters 
# 
P<-1       #par value 

# Define Vasicek Model Parameters 
# 
rev.rate<-0.3     #speed of reversion 
long.term.mean<-0.1    #long term level of the mean 
sigma<-0.05     #volatility 
r0<-0.03       #spot interest rate 
Strike<-0.05 
# Define Simulation Parameters 
# 
T<-50       #time to expiry 
ntime<-1000     #number of timesteps 
yearstep<-ntime/T    #yearstep 
npaths<-1000     #number of paths 
dt<-T/ntime      #timestep 
R <- matrix(0,nrow=ntime+1,ncol=npaths) #matrix of Vasicek interest rate values 
B <- matrix(0,nrow=ntime+1,ncol=npaths) # matrix of Bond Prices 

R[1,]<-r0     #specifies that all paths start at specified spot rate 
B[1,]<-P 

# do loop which generates values to fill matrix R with multiple paths of Interest Rates as they evolve over time. 
# stochastic process based on standard normal distribution 

for (j in 1:npaths) { 
for (i in 1:ntime) { 
    dZ <-rnorm(1,mean=0,sd=1)*sqrt(dt) 
    Rij<-R[i,j] 
    Bij<-B[i,j] 
    dr <-rev.rate*(long.term.mean-Rij)*dt+sigma*dZ 
    R[i+1,j]<-Rij+dr 
    B[i+1,j]<-Bij*exp(-R[i+1,j]*dt) 
    } 
} 

t<-seq(0,T,dt) 
par(mfcol = c(3,3)) 


matplot(t, B[,1:pmin(20,npaths)], type="l", lty=1, main="Price of a Zero Coupon Bond", ylab="Price", xlab = "Time to Expiry") 

Answer 1

你举的例子是不可复制的，所以我创造了我希望的结构类似于你一些假的数据。如果这不是您要查找的内容，请告诉我，我会根据需要进行更新。

# Fake data 
ntime <- 100 
mat1 <- matrix(rnorm(ntime+1, 10, 2), nrow=1, ncol=ntime+1) 
mat2 <- matrix(rnorm(ntime+1, 20, 2), nrow=1, ncol=ntime+1) 
mat3 <- matrix(rnorm(ntime+1, 30, 2), nrow=1, ncol=ntime+1) 

matplot(1:(ntime+1), t(mat1), type="l", lty=1, ylim=c(0, max(c(mat1,mat2,mat3))), 
     main="MC Price of a Zero Coupon Bond", 
     ylab="Price", xlab = "Option Exercise Date") 

# Add lines for mat2 and mat3 
lines(1:101, mat2, col="red") 
lines(1:101, mat3, col="blue") 

UPDATE：这是你想要做什么？

matplot(t, t(average.price.at.each.timestep), type="l", lty=1, 
     main="MC Price of a Zero Coupon Bond", ylab="Price", 
     xlab = "Option Exercise Date") 
matlines(t, t(upper.c.l.at), lty=2, col="red") 
matlines(t, t(lower.c.l.at), lty=2, col="green") 

请参见下图。如果您有多个要绘制的列（如更新示例中绘制20个单独路径的示例），并且想要为所有这些列添加较低和较高的CI（尽管这会使绘图不可读），只需使用矩阵分别对应于average.price.at.each.timestep中的每条路径的上限和下限CI值，并使用matlines将它们添加到您现有的多条路径图中。

Answer 2

这是可行的使用ggplot2和reshape2。你有的结构有点尴尬，你可以通过使用数据框而不是矩阵来改善结构。

#Dummy data 
average.price.at.each.timestep <- rnorm(1000, sd=0.01) 
s.e.at.each.time <- rnorm(1000, sd=0.0005, mean=1) 

#CIs (note you can vectorise this): 
upper.c.l.at <- average.price.at.each.timestep+1.96*s.e.at.each.time 
lower.c.l.at <- average.price.at.each.timestep-1.96*s.e.at.each.time 

#create a data frame: 
prices <- data.frame(time = 1:length(average.price.at.each.timestep), price=average.price.at.each.timestep, upperCI= upper.c.l.at, lowerCI= lower.c.l.at) 

library(reshape2) 
#turn the data frame into time, variable, value triplets 
prices.t <- melt(prices, id.vars=c("time")) 

#plot 
library(ggplot2) 
ggplot(prices.t, aes(time, value, colour=variable)) + geom_line() 

这将产生以下情节：

这可以通过使用代替geom_ribbon有所改善：

ggplot(prices, aes(time, price)) + geom_ribbon(aes(ymin=lowerCI, ymax=upperCI), alpha=0.1) + geom_line() 

将会产生该曲线：

Answer 3

下面是另一个稍有不同的ggplot解决方案，不需要您首先计算置信度限制 - ggplot会为您解决。

# create sample dataset 
set.seed(1) # for reproducible example 
B <- matrix(rnorm(1000,mean=rep(10+1:10/2,each=10)),nc=10) 

library(ggplot2) 
library(reshape2) # for melt(...) 
gg <- melt(data.frame(date=1:nrow(B),B), id="date") 
ggplot(gg, aes(x=date,y=value)) + 
    stat_summary(fun.y = mean, geom="line")+ 
    stat_summary(fun.y = function(y)mean(y)-1.96*sd(y)/sqrt(length(y)), geom="line",linetype="dotted", color="blue")+ 
    stat_summary(fun.y = function(y)mean(y)+1.96*sd(y)/sqrt(length(y)), geom="line",linetype="dotted", color="blue")+ 
    theme_bw() 

stat_summary(...)概括了用于X（日期）的给定值的Y值。因此，在第一次调用中，它计算平均值，第二次计算平均值，第二次计算平均值，第三次计算平均值。

您还可以创建一个CL(...)功能，并称之为：

CL <- function(x,level=0.95,type=c("lower","upper")) { 
    fact <- c(lower=-1,upper=1) 
    mean(x) - fact[type]*qnorm((1-level)/2)*sd(x)/sqrt(length(x)) 
} 
ggplot(gg, aes(x=date,y=value)) + 
    stat_summary(fun.y = mean, geom="line")+ 
    stat_summary(fun.y = CL, type="lower", geom="line",linetype="dotted", color="blue")+ 
    stat_summary(fun.y = CL, type="upper", geom="line",linetype="dotted", color="blue")+ 
    theme_bw() 

这将产生相同上面的一个阴谋。