我有以下代码选择 4 行 iris 1000x,并取每个 4 行样本的平均值:
library(dplyr)
iris<- iris
storage<- list()
counter<- 0
for (i in 1:1000) {
# sample 3 randomly selected transects 100 time
tempsample<- iris[sample(1:nrow(iris), 4, replace=F),]
storage[[i]]=tempsample
counter<- counter+1
print(counter)
}
# Unpack results into dataframe
results<- do.call(rbind, storage)
View(results)
results_2<- as.data.frame(results)
results_2<- results_2 %>% mutate(Aggregate = rep(seq(1,ceiling(nrow(results_2)/4)),each = 4))
# View(results_2)
final_results<- aggregate(results_2[,1:4], list(results_2$Aggregate), mean)
# View(final_results)
我想计算每列相对于其真实总体参数的偏差。例如使用SimDesign
's bias()
:
library(SimDesign)
(bias(final_results[,2:5], parameter=c(5,3,2,1), type='relative'))*100
在此代码中,参数的值是假设的真实弹出。数据框中每一列的值。我想以 100 倍的速度执行此过程,以获取数据框中每个变量的偏差估计分布。但是,我不确定如何将所有这些放入 for 循环中(我认为这是要走的路),所以最终输出是一个数据帧,每个 iris 变量都有 100 行偏差测量值。
对此的任何帮助将不胜感激。
#------------------------------------------
更新
尝试为分层样本而不是随机样本运行相同的代码会给我以下错误: *Error in [.data.table
(setDT(copy(iris)), as.vector(sapply(1:1000, function(X) stratified( iris, : i is invalid type (matrix). 也许将来 2 列矩阵可以返回 DT 的元素列表 * 我认为这可能与 setDT 有关?
这是以下代码的结果:
do.call(rbind,lapply(1:100, function(x) {
bias(
setDT(copy(iris))[as.vector(sapply(1:1000, function(X) stratified(iris,group="Species", size=1)))][
, lapply(.SD, mean), by=rep(c(1:1000),4), .SDcols=c(1:4)][,c(2:5)],
parameter=c(5,3,2,1),
type='relative'
)
}))
我研究了使用以下建议的代码:
get_samples <- function(n, sampsize=4) {
rbindlist(lapply(1:n, function(x) {
splitstackshape::stratified(iris, group="Species",sampsize)[, id:=x] }))[
, lapply(.SD, mean), by=.(Species, id)] }
我想我理解这个函数在做什么(选择 4 行分层的虹膜,按物种取每列的平均值),但我不知道如何将它应用于原来的问题(4 * 1000)* 100 来测试偏见(我对此很陌生,如果我遗漏了一些明显的东西,我深表歉意)。