如何解决如何按组/子集进行留一法交叉验证? 第一部分第二部分编辑编辑 2
这个问题是上一个问题 (Linear Regression prediction in R using Leave One out Approach) 的第二部分。
我正在尝试为每个国家/地区构建模型并使用留一法生成线性回归预测。换句话说,在下面的代码中,在构建模型 1 和模型 2 时,使用的“数据”不应是整个数据集。相反,它应该是数据集(国家/地区)的子集。每个国家/地区的数据都应使用基于该国家/地区特定数据构建的模型进行评估。
下面的代码返回一个错误。我如何修改/修复下面的代码来做到这一点?或者有更好的方法吗?
library(modelr)
install.packages("gapminder")
library(gapminder)
data(gapminder)
#CASE 1
model1 <- lm(lifeExp ~ pop,data = gapminder,subset = country)
model2 <- lm(lifeExp ~ pop + gdpPercap,subset = country)
models <- list(fit_model1 = model1,fit_model2 = model2)
gapminder %>% nest_by(continent,country) %>%
bind_cols(
map(1:nrow(gapminder),function(i) {
map_dfc(models,function(model) {
training <- data[-i,]
fit <- lm(model,data = training)
validation <- data[i,]
predict(fit,newdata = validation)
})
}) %>%
bind_rows()
)
解决方法
最简洁直接的解决方案是嵌套 for
循环方法,其中外循环是两个模型公式,而内循环是我们想要省略的统一体。这也可以使用 outer
来完成,我也会在后面展示。
为了清楚起见,我首先展示了如何在每次迭代中省略一个观察(即一行)(第一部分)。我稍后将展示如何省略一个集群(例如国家/地区)(第二部分)。我还使用了内置的 iris
数据集,它更小,因此更容易处理。它包含一个 "Species"
列,用于对应您数据中的 "countries"
。
第一部分
首先,我们将两个公式放入一个列表中,并按照我们希望它们稍后出现在结果列中的方式命名它们。
FOAE <- list(fit1=Petal.Length ~ Sepal.Length,fit2=Petal.Length ~ Sepal.Length + Petal.Width)
对于循环,我们要初始化一个矩阵 im
,其行对应于我们要省略的行数,列对应于模型公式的数量。
im <- matrix(NA,nrow=nrow(iris),ncol=length(FOAE),dimnames=list(NULL,names(FOAE)))
这看起来像这样:
head(im,n=3)
# fit1 fit2
# [1,] NA NA
# [2,] NA NA
# [3,] NA NA
现在我们如上所述循环公式和行。
for (i in seq(FOAE)) {
for(j in seq(nrow(iris))) {
train <- iris[-j,]
test <- iris[j,]
fit <- lm(FOAE[[i]],data=train)
im[j,i] <- predict(fit,newdata=test)
}
}
im
现在已经被填充了,我们可以将它cbind
到原来的iris
数据集来得到我们的结果res1
。
res1 <- cbind(iris,im)
head(res1)
# Sepal.Length Sepal.Width Petal.Length Petal.Width Species fit1 fit2
# 1 5.1 3.5 1.4 0.2 setosa 2.388501 1.611976
# 2 4.9 3.0 1.4 0.2 setosa 2.014324 1.501389
# 3 4.7 3.2 1.3 0.2 setosa 1.639805 1.392955
# 4 4.6 3.1 1.5 0.2 setosa 1.446175 1.333199
# 5 5.0 3.6 1.4 0.2 setosa 2.201646 1.556620
# 6 5.4 3.9 1.7 0.4 setosa 2.944788 2.127184
为了替代 outer
方法,我们将 for
循环内的代码放入 formula
的 Vectorize
中,以便它可以处理矩阵列(即向量)。
FUN1 <- Vectorize(function(x,y) {
train <- iris[-x,]
test <- iris[x,]
fit <- lm(y,data=train)
predict(fit,newdata=test)
})
现在我们将 FOAE
和行 1:nrow(iris)
随后省略,连同 FUN1
放入 outer()
。这已经为我们提供了结果,我们可以以与上述相同的方式 cbind
到 iris
得到我们的结果 res2
。
o1 <- outer(FOAE,1:nrow(iris),FUN1)
res2 <- cbind(iris,o1)
head(res2)
# Sepal.Length Sepal.Width Petal.Length Petal.Width Species fit1 fit2
# 1 5.1 3.5 1.4 0.2 setosa 2.388501 1.611976
# 2 4.9 3.0 1.4 0.2 setosa 2.014324 1.501389
# 3 4.7 3.2 1.3 0.2 setosa 1.639805 1.392955
# 4 4.6 3.1 1.5 0.2 setosa 1.446175 1.333199
# 5 5.0 3.6 1.4 0.2 setosa 2.201646 1.556620
# 6 5.4 3.9 1.7 0.4 setosa 2.944788 2.127184
## test if results are different is negative
stopifnot(all.equal(res1,res2))
第二部分
我们可能会在遗漏集群(即物种或国家/地区)时采用类似的方法。我在这里展示了 outer
方法。我们想要改变的是,我们现在想要忽略属于特定集群的观察,这里 "Species"
(在您的情况下为 "countries"
),我们将这些 unique
值放入向量中Species.u
。由于值采用 "character"
或 "factor"
格式,我们使用 data[!data$cluster %in% x,]
而不是 data[-x,]
对数据进行子集化。因为 predict
会在集群中产生多个值,但我们希望在各自的集群中具有相同的值,所以我们可能想要使用统计数据,例如每个集群的 mean
预测。我们根据集群使用 rownames
。
FUN2 <- Vectorize(function(x,y) {
train <- iris[!iris$Species %in% x,]
test <- iris[iris$Species %in% x,data=train)
mean(predict(fit,newdata=test))
})
Species.u <- unique(iris$Species)
o2 <- `rownames<-`(outer(Species.u,FOAE,FUN2),Species.u)
这现在给了我们一个小于我们的数据集的矩阵。多亏了 rownames
,我们可以match
预测它们所属的集群。
o2
# fit1 fit2
# setosa 3.609943 2.662609
# versicolor 3.785760 3.909919
# virginica 4.911009 5.976922
res3 <- cbind(iris,o2[match(iris$Species,rownames(o2)),])
head(res3)
# Sepal.Length Sepal.Width Petal.Length Petal.Width Species fit1 fit2
# setosa 5.1 3.5 1.4 0.2 setosa 3.609943 2.662609
# setosa.1 4.9 3.0 1.4 0.2 setosa 3.609943 2.662609
# setosa.2 4.7 3.2 1.3 0.2 setosa 3.609943 2.662609
# setosa.3 4.6 3.1 1.5 0.2 setosa 3.609943 2.662609
# setosa.4 5.0 3.6 1.4 0.2 setosa 3.609943 2.662609
# setosa.5 5.4 3.9 1.7 0.4 setosa 3.609943 2.662609
tail(res3)
# Sepal.Length Sepal.Width Petal.Length Petal.Width Species fit1 fit2
# virginica.44 6.7 3.3 5.7 2.5 virginica 4.911009 5.976922
# virginica.45 6.7 3.0 5.2 2.3 virginica 4.911009 5.976922
# virginica.46 6.3 2.5 5.0 1.9 virginica 4.911009 5.976922
# virginica.47 6.5 3.0 5.2 2.0 virginica 4.911009 5.976922
# virginica.48 6.2 3.4 5.4 2.3 virginica 4.911009 5.976922
# virginica.49 5.9 3.0 5.1 1.8 virginica 4.911009 5.976922
编辑
在这个版本的 FUN2
,FUN3
中,每个集群的模型的输出都经过rbind
ed(当然是两列,因为有两个模型)。
FUN3 <- Vectorize(function(x,data=train)
(predict(fit,newdata=test))
},SIMPLIFY=F)
Species.u <- unique(iris$Species)
o3 <- `rownames<-`(outer(Species.u,FUN3),Species.u)
res32 <- cbind(iris,apply(o3,2,unlist))
head(res32)
# Sepal.Length Sepal.Width Petal.Length Petal.Width Species fit1 fit2
# setosa.1 5.1 3.5 1.4 0.2 setosa 3.706940 2.678255
# setosa.2 4.9 3.0 1.4 0.2 setosa 3.500562 2.547587
# setosa.3 4.7 3.2 1.3 0.2 setosa 3.294183 2.416919
# setosa.4 4.6 3.1 1.5 0.2 setosa 3.190994 2.351586
# setosa.5 5.0 3.6 1.4 0.2 setosa 3.603751 2.612921
# setosa.6 5.4 3.9 1.7 0.4 setosa 4.016508 3.073249
编辑 2
正如我在您的评论中了解到的,您需要 1. 沿着集群的数据子集。这将是下面 ss
中的 FUN4
。然后 ss
也通过在子集 z
的行上省略一行 ss
进行子集化。
FUN4 <- Vectorize(function(x,y) {
## subsets first by cluster then by row
ss <- iris[iris$Species %in% x,] ## cluster subset
sapply(1:nrow(ss),function(z) { ## subset rows using `sapply`
train <- ss[-z,] ## train data w/o row z
test <- ss[z,] ## test data for `predict`,just row z
fit <- lm(y,data=train)
predict(fit,newdata=test)
})
},SIMPLIFY=F)
## the two models
FOAE <- list(fit1=Petal.Length ~ Sepal.Length,fit2=Petal.Length ~ Sepal.Length + Petal.Width)
## unique cluster names
Species.u <- unique(iris$Species)
## with the `outer` we iterate over all the permutations of clusters and models `FOAE`.
o4 <- `rownames<-`(outer(Species.u,FUN4),Species.u)
## `unlist`ed result is directly `cbind`able to original data
res4 <- cbind(iris,apply(o4,unlist))
## result
head(res4)
# Sepal.Length Sepal.Width Petal.Length Petal.Width Species fit1 fit2
# setosa.1 5.1 3.5 1.4 0.2 setosa 1.476004 1.451029
# setosa.2 4.9 3.0 1.4 0.2 setosa 1.449120 1.431737
# setosa.3 4.7 3.2 1.3 0.2 setosa 1.426185 1.416492
# setosa.4 4.6 3.1 1.5 0.2 setosa 1.404040 1.398103
# setosa.5 5.0 3.6 1.4 0.2 setosa 1.462460 1.441295
# setosa.6 5.4 3.9 1.7 0.4 setosa 1.504990 1.559045
版权声明:本文内容由互联网用户自发贡献,该文观点与技术仅代表作者本人。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如发现本站有涉嫌侵权/违法违规的内容, 请发送邮件至 dio@foxmail.com 举报,一经查实,本站将立刻删除。