我有以下代碼選擇 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 來測驗偏見（我對此很陌生，如果我遺漏了一些明顯的東西，我深表歉意）。
uj5u.com熱心網友回復：
這是一種方法。我對您的代碼進行了一些小的更改，并將其包裝在一個函式中。然后，使用lapply一個序列說1:10or 1:100，每次運行你的函式，并將結果bias從SimDesign包中提供給你的函式。然后行系結結果串列
library(dplyr)
get_samples <- function(df, size=4, n=1000) {
storage<- list()
counter<- 0
for (i in 1:1000) {
tempsample<- df[sample(1:nrow(df), size, replace=F),]
storage[[i]]=tempsample
counter<- counter 1
}
results<- do.call(rbind, storage)
results_2<- as.data.frame(results)
results_2<- results_2 %>% mutate(Aggregate = rep(seq(1,ceiling(nrow(results_2)/size)),each = size))
final_results<- aggregate(results_2[,1:size], list(results_2$Aggregate), mean)
return(final_results)
}
iris=iris
replicates = lapply(1:10, function(x) {
result = get_samples(iris)
(bias(result[,2:5], parameter=c(5,3,2,1), type='relative'))*100
})
replicates = do.call(rbind, replicates)
輸出：
Sepal.Length Sepal.Width Petal.Length Petal.Width
[1,] 41.50617 3.292500 86.77408 8.859333
[2,] 43.26058 2.763500 90.20758 10.825917
[3,] 43.46642 3.551750 90.11767 10.576250
[4,] 41.94683 2.970833 86.89625 8.817000
[5,] 42.08733 3.380917 86.78642 8.996667
[6,] 42.13050 2.942250 88.02983 9.707500
[7,] 43.07383 2.775500 89.04583 10.102083
[8,] 44.10192 2.895167 91.27208 11.188500
[9,] 41.29408 2.314750 87.59208 9.244333
[10,] 42.77450 2.781583 90.37342 10.789500
快速解決問題
library(SimDesign)
library(data.table)
do.call(rbind,lapply(1:100, function(x) {
bias(
setDT(copy(iris))[as.vector(sapply(1:1000, function(X) sample(1:nrow(iris),4)))][
, lapply(.SD, mean), by=rep(c(1:1000),4), .SDcols=c(1:4)][,c(2:5)],
parameter=c(5,3,2,1),
type='relative'
)
}))
uj5u.com熱心網友回復：
由于您正在使用mutate，您可以考慮使用tidyverse.
map_df(1:1000, ~ sample_n(iris, 4, replace = FALSE)) %>%
glimpse() %>%
mutate(Aggregate_col = rep(seq(1, ceiling(n() / 4)), each = 4)) %>%
glimpse() %>%
select(starts_with("Sepal"),
starts_with("Petal"),
matches("Aggregate")) %>%
group_by(Aggregate_col) %>%
summarise(across(.cols = everything(), ~ mean(.x, na.rm = TRUE)))
筆記：
在下面的示例中，您的第一個回圈被替換為：
map_df(1:1000, ~ sample_n(iris, 4, replace = FALSE))map_x可以用來迭代一個串列，或者在這種情況下是一個整數向量1:1000，如果唯一的目的是重復呼叫函式，并將結果系結到所需的格式，在這種情況下是 adata.frame。您可以
glimpse在資料轉換管道中利用 while 來避免View重復呼叫select提供一種按名稱或部分匹配選擇列的可讀方式。這通常比在添加/洗掉變數時按索引選擇列更安全
轉載請註明出處，本文鏈接：https://www.uj5u.com/yidong/433026.html