2017-03-16 93 views
1

我试图绘制对同一地块不同的方程一些累积分布函数,这是我写的代码:情节没有显示CDF

library (triangle) 
library(lattice) 
library(latticeExtra) 
n = 1000 

W1 = rtriangle(n,3128,3250) 
W2 = rtriangle(n,3065,3149) 
SO = rtriangle(n,0.2,0.3) 

MCtab <- data.frame(W1,W2,SO) 

set.seed(1) 
for (n in 1:n) { 
NPV30 <- (1*W1 + 2*W2 + 3*SO)} 

set.seed(1) 
for (n in 1:n) { 
NPV50 <- ecdf((4*W1 + 5*W2 + 6*SO))} 
set.seed(1) 
for (n in 1:n) { 
NPV70 <- ecdf((7*W1 + 8*W2 + 9*SO))} 

MCtab2 <- data.frame(NPV30,NPV50,NPV70) 



plot(NPV30, verticals=TRUE, main= 'Polymer NPV CDF', do.points=FALSE, col='red') 
plot(NPV50, verticals=TRUE, do.points=FALSE, add=TRUE, col='brown') 
plot(NPV70, verticals=TRUE, do.points=FALSE, add=TRUE, col='orange') 

只有一个曲线如图不能似乎找出原因,或者如果有人有更好的方法会很好。谢谢

+1

一种可能性是,后来的值太大初始待观察情节维度。您可以添加包含所有数据点的最小值和最大值的ylim参数。 – lmo

+0

@lmo谢谢,这似乎是最有可能的原因 – John

+0

对我的答案的任何反馈? –

回答

0

我假设在示例代码中NPV30赋值之前缺少ecdf。我认为它没有显示值的实际原因是plot.ecdf方法根本不支持add参数。但你可以很容易地推出你自己的。

您为这些函数选择的值给出了不同的分布,累积分布图很无聊 - 只是矩形。所以我调整了NPV50NPV70的定义值,使三个分布更加相似。

因此,这里是我使用的代码:

library (triangle) 
library(lattice) 
library(latticeExtra) 
n = 50 

W1 = rtriangle(n,3128,3250) 
W2 = rtriangle(n,3065,3149) 
SO = rtriangle(n,0.2,0.3) 

MCtab <- data.frame(W1,W2,SO) 

set.seed(1) 
for (n in 1:n) { 
    NPV30 <- ecdf((1*W1 + 2*W2 + 3*SO))} 

set.seed(1) 
for (n in 1:n) { 
    NPV50 <- ecdf((1*W1 + 2.2*W2 + 3.1*SO))} 
set.seed(1) 
for (n in 1:n) { 
    NPV70 <- ecdf((1*W1 + 2.4*W2 + 3.2*SO))} 

MCtab2 <- data.frame(NPV30,NPV50,NPV70) 

qNPV30 <- quantile(NPV30) 
qNPV50 <- quantile(NPV50) 
qNPV70 <- quantile(NPV70) 
xmin <- min(c(qNPV30[1],qNPV50[1],qNPV70[1])) 
xmax <- max(c(qNPV30[5],qNPV50[5],qNPV70[5])) 
xdlt <- xmax-xmin 

npts <- 200 
x <- (1:npts)*xdlt/npts + xmin 
y30 <- NPV30(x) 
y50 <- NPV50(x) 
y70 <- NPV70(x) 

plot(x,y30,main= 'Polymer NPV CDF',xlim=c(xmin,xmax),type='l',col='red') 
lines(x,y50,xlim=c(xmin,xmax),type='l',col='blue',add=T) 
lines(x,y70,xlim=c(xmin,xmax),type='l',col='green',add=T) 

#plot(NPV30, verticals=TRUE, main= 'Polymer NPV CDF', do.points=FALSE, col='red') 
#plot(NPV50, verticals=TRUE, do.points=F, add=TRUE, col='blue') 
#plot(NPV70, verticals=TRUE, do.points=F, add=TRUE, col='green') 

,这里是输出结果:

enter image description here

+0

感谢您的澄清(道歉为非常迟到的答复) – John