saving output from a nested for-loop in R

Question

I'm trying to figure out how to store the output of a simulation that I will run for many iterations. I have a simplified example here where I'm simulating changes in winter temperature over a 50-year period, and I'd like to do this 3 times:

# function to draw the winter temperature in a particular year
draw_winter_temperature <- function(){
  intercept <- winterIntercept
  beta1 <- rnorm(1, winterBeta, winterBetaSD)
  val = rnorm(1, intercept + beta1 * years[i], winterSigma)
  return(val)
}

# parameters for winter temperature function
winterIntercept <- -0.633
winterBeta <- 0.08
winterBetaSD <- 0.04
winterSigma <- 1.30 

# number of iterations
iter = seq(1,3,1)

# number of years to run in each simulation
years <- seq(1,50,1)

# matrix to store 50 winter temperatures for each iteration? Not sure best approach
temps <- matrix(NA, nrow = length(years), ncol = length(iter))

# start of nested for loop - loop over each iteration
for (i in seq_along(iter)) {

  # within each iteration, loop over each year
  for (j in seq_along(years)) {

  # draw winter temp. anomaly
  winterTemp <- draw_winter_temperature()

  # store winter temps, not sure what to do here
  temps <- winterTemp

  }
}

I guess I'm a little unsure about the indexing here to save the output in my empty matrix (or to save it in another format). I should say that there are other functions and data being stored within my loop over years, so I'm generally confused about the best way to store data from each of the iterations.

neilfws · Accepted Answer · 2019-03-04 03:53:27Z

I'd avoid the loops altogether. Create a data frame with the iterations and years:

data1 <- data.frame(iter = rep(1:3, each = 50), 
                    year = rep(1:50, 3))

Define the function a little differently, so as all the default parameters are supplied:

draw_winter_temperature <- function(intercept = -0.633, 
                                    winterbeta = 0.08, 
                                    winterbetasd = 0.04,
                                    wintersigma = 1.30,
                                    year = 1){
  beta1 <- rnorm(1, winterbeta, winterbetasd)
  val   <- rnorm(1, intercept + beta1 * year, wintersigma)
  val
}

Use purrr::map_dbl() and dplyr::mutate() to create a column with sampled values for each year:

library(dplyr)
library(purrr)

data1 <- data1 %>% 
  mutate(winterTemp  = map_dbl(year, ~draw_winter_temperature(year = .)))

Santiago Capobianco · Accepted Answer · 2019-03-04 19:07:36Z

I think you need the [j, i] subset in temps. Also, redefine the function to move across the years argument:

> rm(list = ls())
> 
> # parameters for winter temperature function
> winterIntercept <- -0.633
> winterBeta <- 0.08
> winterBetaSD <- 0.04
> winterSigma <- 1.30 
> 
> # number of iterations
> iter = seq(1,3,1)
> 
> # number of years to run in each simulation
> years <- seq(1,50,1)
> 
> # matrix to store 50 winter temperatures for each iteration? Not sure best approach
> temps <- matrix(NA, nrow = length(years), ncol = length(iter))
> 
> draw_winter_temperature <- function(i) {
+   intercept <- winterIntercept
+   beta1 <- rnorm(1, winterBeta, winterBetaSD)
+   val = rnorm(1, intercept + beta1 * years[i], winterSigma)
+   return(val)
+ }
> 
> 
> # start of nested for loop - loop over each iteration
> for (i in seq_along(iter)) {
+   
+   # within each iteration, loop over each year
+   for (j in seq_along(years)) {
+     
+     # draw winter temp. anomaly
+     winterTemp <- draw_winter_temperature(j)
+     
+     # store winter temps, not sure what to do here
+     temps[j, i] <- winterTemp
+     
+   }
+ }
> 
> temps
             [,1]        [,2]        [,3]
 [1,]  0.07484416 -2.09659305  1.36146643
 [2,]  0.91077388  0.91630545  0.15919461
 [3,]  2.07646459  1.47967364 -0.80436293
 [4,] -1.06727442 -1.35760417 -2.23396119
 [5,]  0.67833459 -0.49733006 -0.08458830
 [6,] -2.94754201 -0.83298461 -0.09547269
 [7,] -0.27718212 -1.26709144  1.07166328
 [8,] -2.48913267  1.15023058 -0.45768205
 [9,]  0.81201831  1.73791432 -1.12876304
[10,] -1.26716059 -0.59580422 -1.28796699
[11,]  1.73244202  1.45400888 -0.13374732
[12,]  0.06181381  1.42056333 -0.28493854
[13,]  1.32197105 -0.26822322 -0.64071024
[14,]  0.31681581  2.17959219  1.05139421
[15,] -0.40649165  0.40584711 -1.24801466
[16,] -0.37852770 -0.01912727 -1.57738425
[17,]  1.07342597 -0.37121912 -1.78913387
[18,]  2.37652192  2.34135284  1.57626300
[19,]  2.26526300 -0.22457651  1.22608845
[20,]  0.43476700  1.20513841  3.15556064
[21,]  2.59591314 -0.30772162 -0.18688845
[22,]  1.06465892  1.71525510 -1.96024033
[23,]  0.33484399  1.81225896  2.11832405
[24,]  1.90063594  3.51556209  1.07306117
[25,] -0.01940422  2.28460208  3.21963212
[26,]  3.28816991  2.95440512  3.24107428
[27,]  0.96753190 -0.10513793 -1.21244216
[28,]  2.54233795 -0.71844973  1.81812858
[29,]  2.17954621  2.98536138  3.66748584
[30,]  0.81877427 -1.35965678  2.56338093
[31,]  0.29897273  1.11832993  0.87086081
[32,]  3.41605791  4.72562882  4.07046433
[33,]  3.48441756 -0.09328612  1.64075691
[34,] -1.45902220 -0.38229884  4.62577573
[35,]  0.76568048  3.05618914 -2.04568644
[36,]  0.51939827  0.46442133  2.45079102
[37,]  3.31122478  3.72294771  3.79475270
[38,]  3.22139644  2.33785724 -0.38787472
[39,]  2.95662591  1.14475461  2.97768620
[40,]  0.41263989  3.12365715 -0.05141565
[41,]  1.38138984  3.76526281  1.96114126
[42,] -0.18932163  3.05150559  1.01252868
[43,]  2.60596068  4.58785954 -0.62515705
[44,]  4.11088547  1.61733269  2.51478319
[45,]  2.54634136  5.03258017  3.69739035
[46,]  3.90668226  5.10453707  9.44805817
[47,] -0.01088235 -0.72829163  5.23958084
[48,]  5.69817388  0.12111318  4.43908874
[49,]  0.62863949  1.44337538  6.37749772
[50,]  6.03095905 -1.06661804  2.96625244

It can also be done with nested sapply:

> sapply(seq(1, 3, 1), function(x){
+   sapply(seq(1, 50, 1), function(i) {draw_winter_temperature(i)})
+ })
             [,1]       [,2]        [,3]
 [1,] -0.89724777 -0.3248484 -0.90899374
 [2,] -1.94179327 -1.1282839  0.31149541
 [3,] -0.42025521  0.7724711  0.74888190
 [4,] -1.33109020  1.9129040  0.87848969
 [5,] -0.47958978  0.0847167  0.41101536
 [6,] -0.79725942 -0.3091579 -1.24232376
 [7,]  1.77518276 -0.1922316 -1.47421230
 [8,]  0.69493304 -0.2019549  0.13570145
 [9,] -1.30457205  0.3104541  0.42987422
[10,]  0.17036195 -1.7092747 -0.83657188
[11,]  0.32147101 -0.9007392 -1.19378838
[12,]  0.67253811  1.3546455  0.65988074
[13,]  1.50350007  0.5492711  1.79631547
[14,] -2.51666246 -0.9415391 -0.29810167
[15,] -0.40993616 -0.4815619 -0.64542730
[16,]  1.77118058  0.1315932 -0.75448968
[17,]  0.50048465  2.7956318 -0.87173764
[18,] -3.21569892 -0.3323585 -1.23359412
[19,]  0.31825702  0.3584551  0.06093396
[20,]  3.89144244 -2.4149095  1.78682882
[21,] -0.43390072  0.7793848 -0.10369739
[22,] -2.86730666 -0.7670008  1.52314298
[23,] -1.79228394  0.9572841  0.48042120
[24,] -0.43718475 -0.1275890  1.21647106
[25,] -1.15994958  1.6406028  4.60036997
[26,]  0.01317481  1.4994584  2.66739115
[27,]  4.17871392  0.7907429  2.73681377
[28,]  4.28581141  1.6371291  0.68691671
[29,]  1.93210763  2.9950647  4.21374229
[30,]  2.58606149  1.0610567 -0.85244266
[31,]  1.80716913  4.5783347  2.52572309
[32,] -0.59753770  1.4036211  4.67023114
[33,] -0.62456359  3.8556583  1.45738330
[34,]  3.97709081  1.9414286  2.13243487
[35,] -0.44078059  0.1845302  2.93266443
[36,]  3.95186775  0.6522162  1.33296220
[37,]  2.44364830  2.8547030  1.20660793
[38,]  4.41777409  1.9200504  2.71936027
[39,]  0.28099830  2.7876076  3.92335543
[40,]  3.40916089  0.1663190  0.42249808
[41,]  4.14691320  2.2004945  2.03032191
[42,]  0.61603992  4.6339235  2.82358688
[43,]  3.15565113  4.2119661  5.04334175
[44,]  2.63957253  5.9546210  2.77191218
[45,] -1.50822755  0.9417640  6.71999760
[46,]  1.29101617  3.2551090  1.59456129
[47,]  5.98411236  2.0040345  4.03826529
[48,]  2.33718663 -3.1029647  2.84654953
[49,] -1.77860722  4.4929658  2.93310143
[50,]  7.12512178  5.3284199  2.97240287

Hope it helps.

@ Santiago Capobianco I tried your first recommendation, but in that case the temperatures do not increase over time as they should. I should've also mentioned that there are many other functions and things happening inside the for-loop, so I'm not sure if sapply would work for my situtation.
I've just realize that the function must be redefined like in the second example. I'll change it in a sec.

Collectives™ on Stack Overflow

saving output from a nested for-loop in R

2 Answers 2

Comments

2 Comments

Your Answer

Hot Network Questions

Collectives™ on Stack Overflow

2 Answers 2

Comments

2 Comments

Your Answer

Sign up or log in

Post as a guest

Related