Problems with appending t.test results in a for loop

user_8114807 0
Viewed 44

Let me take simulated datasets to explain:

I have dataset dt and dt1

# dataset 1 `dt`
set.seed(12)
dt <- rnorm(5000,mean=10,sd=1)
dt <- data.frame(dt)
dt$group <- c("case","control")
colnames(dt) <- c("severity","group")
head(dt)

severity   group
1  8.519432    case
2 11.577169 control
3  9.043256    case
4  9.079995 control
5  8.002358    case
6  9.727704 control


# dataset 2 `dt2`
set.seed(12)
dt2 <- rnorm(200,mean=12,sd=1)
dt2 <- data.frame(dt2)
dt2$group <- c("case2","control2")
colnames(dt2) <- c("severity","group")
head(dt2)

  severity    group
1 10.51943    case2
2 13.57717 control2
3 11.04326    case2
4 11.07999 control2
5 10.00236    case2
6 11.72770 control2

I am building one 1000 iterations for loop to do the following steps:

  • randomly take 500 rows from the dt and save as dt_sub
  • rbind dt_sub with dt2 and save as bd
  • select only rows with group as either case2 or control from the bd dataset (only cares the difference between these two groups)
  • t.tests on the variable severity between the case2 and control group
  • output t.tests results to t
  • use a for loop to repeat 1000 times
  • iteratively appends all t.test results to a dataframe results

Following is the code that I built in r

library(broom)
library(dplyr)

iter <- 1000
t <- data.frame()
for (i in 1:iter) {
  dt_sub <- dt[sample(nrow(dt),500),]
  bd <- rbind(dt_sub,dt2)
  compare <- filter(bd, group %in% c("case2", "control"))
  compare %>% group_by(group) %>% do(tidy(t.test(severity ~ group,data = compare))) -> t
  t$iter <- i
}

  results <- do.call(rbind,t)

My question is, this code works well when iter=1, but how should I set the compare %>% group_by(group) %>% do(tidy(t.test(severity ~ group,data = compare))) -> t line to ensure each run's t.test results will not be overwritten when iter ≥ 1? I tried t[i] but failed, anyone could advise please?

Thanks.

append iteration r for-loop
1 Answers

Create a function which runs the process once.

library(broom)
library(dplyr)

t_test_function <- function() {
  dt_sub <- dt[sample(nrow(dt),500),]
  bd <- rbind(dt_sub,dt2)
  compare <- filter(bd, group %in% c("case2", "control"))
  compare %>% 
    group_by(group) %>% 
    do(tidy(t.test(severity ~ group,data = compare))) %>%
    ungroup
}
t_test_function()

#  group estimate estimate1 estimate2 statistic  p.value parameter conf.low
#  <chr>    <dbl>     <dbl>     <dbl>     <dbl>    <dbl>     <dbl>    <dbl>
#1 case2     1.94      11.9      9.99      17.4 9.40e-42      199.     #1.72
#2 cont…     1.94      11.9      9.99      17.4 9.40e-42      199.     1.72
# … with 3 more variables: conf.high <dbl>, method <chr>,
#   alternative <chr>

Now you can call this iter times using replicate and combine the dataset.

iter <- 5
results <- bind_rows(replicate(iter, t_test_function(), simplify = FALSE), .id = 'iter')

# A tibble: 10 x 12
#   iter  group estimate estimate1 estimate2 statistic  p.value parameter
#   <chr> <chr>    <dbl>     <dbl>     <dbl>     <dbl>    <dbl>     <dbl>
# 1 1     case2     1.88      11.9     10.1       17.3 1.05e-40      189.
# 2 1     cont…     1.88      11.9     10.1       17.3 1.05e-40      189.
# 3 2     case2     1.96      11.9      9.97      17.8 9.88e-43      194.
# 4 2     cont…     1.96      11.9      9.97      17.8 9.88e-43      194.
# 5 3     case2     1.94      11.9      9.99      17.9 3.76e-42      184.
# 6 3     cont…     1.94      11.9      9.99      17.9 3.76e-42      184.
# 7 4     case2     2.03      11.9      9.90      18.6 1.82e-44      189.
# 8 4     cont…     2.03      11.9      9.90      18.6 1.82e-44      189.
# 9 5     case2     1.96      11.9      9.97      18.1 7.05e-43      187.
#10 5     cont…     1.96      11.9      9.97      18.1 7.05e-43      187.
# … with 4 more variables: conf.low <dbl>, conf.high <dbl>, method <chr>,
#   alternative <chr>
user_39629140
Related