1

本质上,我想建立一个用于敏感性分析的蜘蛛图。我想将我的数据分成 10 个部分,并找到每个部分的平均结果值(在第 4 列中)。应根据每个变量列中数据的第 10、第 20、第 30、第 40 等百分位数选择批次。我得到了这个工作,但我认为必须有一个更简单的方法来做到这一点。

我的代码:

##Make some data and put it into a matrix.

c <- 1000
v1 <- rnorm (c, 100, 15)
v2 <- rnorm (c, 80, 10)
v3 <- rnorm (c, 50, 5)
r1 <- ((v1*v2^2)/v3)
data <- cbind (v1,v2)
data <- cbind (data, v3)
data <- cbind (data, r1)

##Sort matrix by first column.
data <- as.matrix(data[order(data[,1]),])

##Find mean of column 4 values corresponding to the smallest 10% (and 20%, and 30%,     etc.) of column 1 values.
a1 <- mean (data[1:(c/10),4])
a2 <- mean (data[(c/10):(2*c/10),4])
a3 <- mean (data[(2*c/10):(3*c/10),4])
a4 <- mean (data[(3*c/10):(4*c/10),4])
a5 <- mean (data[(4*c/10):(5*c/10),4])
a6 <- mean (data[(5*c/10):(6*c/10),4])
a7 <- mean (data[(6*c/10):(7*c/10),4])
a8 <- mean (data[(7*c/10):(8*c/10),4])
a9 <- mean (data[(8*c/10):(9*c/10),4])
a10 <- mean (data[(9*c/10):c,4])

##Combine into a vector.
a <- as.vector(c(a1, a2, a3, a4, a5, a6, a7, a8, a9, a10))

##Repeat for data sorted by columns 2 and 3 respectively.
data <- as.matrix(data[order(data[,2]),])

a1 <- mean (data[1:(c/10),4])
a2 <- mean (data[(c/10):(2*c/10),4])
a3 <- mean (data[(2*c/10):(3*c/10),4])
a4 <- mean (data[(3*c/10):(4*c/10),4])
a5 <- mean (data[(4*c/10):(5*c/10),4])
a6 <- mean (data[(5*c/10):(6*c/10),4])
a7 <- mean (data[(6*c/10):(7*c/10),4])
a8 <- mean (data[(7*c/10):(8*c/10),4])
a9 <- mean (data[(8*c/10):(9*c/10),4])
a10 <- mean (data[(9*c/10):c,4])

b <- as.vector(c(a1, a2, a3, a4, a5, a6, a7, a8, a9, a10))

data <- as.matrix(data[order(data[,3]),])

a1 <- mean (data[1:(c/10),4])
a2 <- mean (data[(c/10):(2*c/10),4])
a3 <- mean (data[(2*c/10):(3*c/10),4])
a4 <- mean (data[(3*c/10):(4*c/10),4])
a5 <- mean (data[(4*c/10):(5*c/10),4])
a6 <- mean (data[(5*c/10):(6*c/10),4])
a7 <- mean (data[(6*c/10):(7*c/10),4])
a8 <- mean (data[(7*c/10):(8*c/10),4])
a9 <- mean (data[(8*c/10):(9*c/10),4])
a10 <- mean (data[(9*c/10):c,4])

d <- as.vector(c(a1, a2, a3, a4, a5, a6, a7, a8, a9, a10))

##Make a pretty chart
plot (a, type = "o", col = "red") 
lines (b, type = "o", col = "blue") 
lines (d, type = "o", col = "green") 
4

2 回答 2

2

这是一些做同样事情的代码,但对 R 来说更紧凑和惯用。

n <- 1000
# changed from c to n since you use c again later as something else
v1 <- rnorm (n, 100, 15)
v2 <- rnorm (n, 80, 10)
v3 <- rnorm (n, 50, 5)
r1 <- ((v1*v2^2)/v3)

DF <- data.frame(v1, v2, v3, r1)
# A data.frame seems like it would be a better fit for this

library("Hmisc")
# The Hmisc package has a function which splits in to quantiles, so use it
DF <- transform(DF, 
                v1.decile = cut2(v1, g=10),
                v2.decile = cut2(v2, g=10),
                v3.decile = cut2(v3, g=10))
# add three new variables to the data frame which indicate which decile each
# value belongs to, for each of v1, v2, and v3
a <- aggregate(DF$r1, list(DF$v1.decile), mean)$x
# why add the new variables? because aggregate can perform an operation on
# groups of one variable defined by the value of another variable
b <- aggregate(DF$r1, list(DF$v2.decile), mean)$x
c <- aggregate(DF$r1, list(DF$v3.decile), mean)$x

然后你可以像以前一样制作情节。


编辑:

Ananda Mahto 的回答指出了我忘记的聚合函数的函数版本。aggregate你可以更清楚地写出这些线条

a <- aggregate(r1 ~ v1.decile, DF, mean)$r1
b <- aggregate(r1 ~ v2.decile, DF, mean)$r1
c <- aggregate(r1 ~ v3.decile, DF, mean)$r1
于 2013-01-30T05:55:40.467 回答
1

这在概念上与 Brian Diggs 的答案非常相似,但不依赖于您的输入是 adata.frame也不依赖于加载任何包。它还介绍了matplot,它将为您提供绘图,而不必一次绘制每一列。

这是您的数据:

set.seed(1) # make it reproducible 
n <- 1000
v1 <- rnorm (c, 100, 15)
v2 <- rnorm (c, 80, 10)
v3 <- rnorm (c, 50, 5)
r1 <- ((v1*v2^2)/v3)
data <- cbind (v1, v2, v3, r1)
rm(v1, v2, v3, r1) # Cleanup

head(data)
#             v1       v2       v3        r1
# [1,]  90.60319 95.11781 54.59489 15014.651
# [2,] 102.75465 83.89843 53.91068 13416.349
# [3,]  87.46557 73.78759 50.37282  9453.824
# [4,] 123.92921 57.85300 40.05324 10355.899
# [5,] 104.94262 91.24931 53.09913 16455.977
# [6,]  87.69297 79.55066 49.71936 11161.612

我们将使用它sapply来执行我们的聚合。这将产生一个我们可以轻松绘制的矩阵。

myAggVars <- c("v1", "v2", "v3")
temp <- sapply(myAggVars, function(x) {
  aggregate(r1 ~ cut(get(x), quantile(get(x), probs = seq(0, 1, .1)), 
                     include.lowest = TRUE), data, mean)[[2]]
})
temp
#              v1        v2        v3
#  [1,]  9453.824 10355.899 10355.899
#  [2,] 11161.612  9453.824 20834.485
#  [3,] 15014.651 11161.612 17755.902
#  [4,] 13528.961 13896.830 13896.830
#  [5,] 13416.349 13416.349 11161.612
#  [6,] 16455.977 13528.961  9453.824
#  [7,] 13896.830 17755.902 13528.961
#  [8,] 17755.902 20834.485 16455.977
#  [9,] 20834.485 16455.977 13416.349
# [10,] 10355.899 15014.651 15014.651

这是绘图步骤:

matplot(temp, type = "o", pch = 1)

结果:

在此处输入图像描述

于 2013-01-30T07:41:24.577 回答