[R] Monte Carlo on simple regression

[Eric Berger]

Hi Ogbos,
I took a closer look at your code.
Here's a modified version (using dummy data) that seems to do what you want.
Hopefully this will make it clear what you need to to.

nn <- 100
lDf <- data.frame(Li=rnorm(nn),CR=rnorm(nn))

fit<-lm(Li~CR, data=lDf)
a<-summary(fit)

N <- nrow(lDf)
C <- 50         # desired number of subsamples
S <- 38         # desired sample size

sumb2 <- 0
for (i in 1:C){   # a loop over the number of subsamples
  set.seed(3*i)   # a different seed for each subsample
  subsample <- lDf[sample(1:N, size=S, replace=TRUE), ]
  mod <- lm(Li~CR,data=subsample)
  #sum b2 for all subsamples:
  sumb2 <- sumb2 + coef(mod)[[2]]
}
print(sumb2/C, digits = 3)

Best,
Eric



On Wed, Aug 22, 2018 at 7:28 PM, Ogbos Okike <giftedlife2014 using gmail.com>
wrote:

> Hello Erick,
>
> Thanks again.
> Another line indicated error:
>
> source("script.R")
> Error in eval(predvars, data, env) :
>   numeric 'envir' arg not of length one
> Thank you for additional assitance.
> Ogbos
>
>
>
> On Wed, Aug 22, 2018 at 5:23 PM Eric Berger <ericjberger using gmail.com> wrote:
>
>> You have an extra comma ... it should be
>>
>> Li[sample(1:N, size = S, replace = TRUE)]
>>
>> i.e. no comma after the closing parenthesis
>>
>>
>>
>> On Wed, Aug 22, 2018 at 7:20 PM, Ogbos Okike <giftedlife2014 using gmail.com>
>> wrote:
>>
>>> Hello Eric,
>>> Thanks for this.
>>>
>>> I tried it. It went but another problem prevents the code from running.
>>>  source("script.R")
>>> Error in Li[sample(1:N, size = S, replace = TRUE), ] :
>>>   incorrect number of dimensions
>>>
>>> The error is coming from the line:
>>>  subsample <- Li[sample(1:N, size=S, replace=TRUE), ]
>>>
>>> I tried to replace Li with N but it didn't go. I also tried replacing it
>>> with length(Li). The same error remains.
>>>
>>> Thank so much for looking at this again.
>>>
>>> Ogbos
>>>
>>>
>>> On Wed, Aug 22, 2018 at 5:06 PM Eric Berger <ericjberger using gmail.com>
>>> wrote:
>>>
>>>> Li is defined as d1$a which is a vector. You should use
>>>>
>>>> N <- length(Li)
>>>>
>>>> HTH,
>>>> Eric
>>>>
>>>>
>>>> On Wed, Aug 22, 2018 at 6:02 PM, Ogbos Okike <giftedlife2014 using gmail.com>
>>>> wrote:
>>>>
>>>>> Kind R-users,
>>>>> I run a simple regression. I am interested in using the Monte Carlo to
>>>>> test
>>>>> the slope parameter.
>>>>> Here is what I have done:
>>>>> d1<-read.table("Lightcor",col.names=c("a"))
>>>>> d2<-read.table("CRcor",col.names=c("a"))
>>>>>  Li<-d1$a
>>>>> CR<-d2$a
>>>>>
>>>>>  fit<-lm(Li~CR)
>>>>>  a<-summary(fit)
>>>>> a gives the slope as 88.15
>>>>>
>>>>> Problem: I now what to repeat the samples to access this coefficient.
>>>>> Following one of the related examples I got online, I did (tried to
>>>>> modify):
>>>>>
>>>>> N <- nrow(Li) # returns the number of observations in the dataset
>>>>> C <- 50         # desired number of subsamples
>>>>> S <- 38         # desired sample size
>>>>>
>>>>> sumb2 <- 0
>>>>> for (i in 1:C){   # a loop over the number of subsamples
>>>>>   set.seed(3*i)   # a different seed for each subsample
>>>>>   subsample <- Li[sample(1:N, size=S, replace=TRUE), ]
>>>>>   mod <- lm(Li~CR,data=subsample)
>>>>>   #sum b2 for all subsamples:
>>>>>   sumb2 <- sumb2 + coef(mod)[[2]]
>>>>> }
>>>>> print(sumb2/C, digits = 3)
>>>>>
>>>>>    But when I run the script, I had error message:
>>>>> Error in 1:N : argument of length 0
>>>>> My data:
>>>>> Li        CR
>>>>> 74281 8449
>>>>> 92473 8148
>>>>> 62310 8520
>>>>> 71219 8264
>>>>> 33469 8389
>>>>> 75768 7499
>>>>> 61636 7821
>>>>> 103829 8468
>>>>> 87336 8568
>>>>> 129443 8190
>>>>> 97682 8539
>>>>> 106918 8502
>>>>> 97171 8578
>>>>> 48012 8181
>>>>> 93086 8631
>>>>> 92374 8562
>>>>> 113010 8404
>>>>> 66956 8592
>>>>> 133037 8632
>>>>> 108849 8644
>>>>> 81544 8442
>>>>> 105072 8615
>>>>> 143437 7724
>>>>> 153294 7829
>>>>> 123735 8682
>>>>> 154738 8756
>>>>> 100760 8839
>>>>> 108034 8839
>>>>> 81826 8858
>>>>> 116901 8847
>>>>> 80780 8869
>>>>> 122684 8736
>>>>> 141716 9087
>>>>> 144315 9166
>>>>> 162078 9147
>>>>> 163184 9267
>>>>> 150688 9275
>>>>> 200848 9259
>>>>> 221570 8943
>>>>> 192424 8564
>>>>> 173024 9282
>>>>> 197326 9318
>>>>> 209344 9293
>>>>> 220201 9242
>>>>> 212626 9324
>>>>> 218115 9319
>>>>> 170001 9314
>>>>> 187490 9346
>>>>> 172440 9350
>>>>> 180330 9349
>>>>> 200807 9355
>>>>> 234994 9350
>>>>> 139053 9284
>>>>> 150048 9361
>>>>> 203650 9346
>>>>> 233331 9369
>>>>> 198790 9340
>>>>> 164060 9382
>>>>> 198000 9401
>>>>> 201707 9355
>>>>> 179257 9369
>>>>> 188736 9298
>>>>> 243392 9393
>>>>> 246040 9374
>>>>> 269058 9364
>>>>> 201657 9370
>>>>> 187942 9354
>>>>> 228514 9305
>>>>> 234000 9392
>>>>> 224431 9395
>>>>> 163502 9398
>>>>> I would be most glad for your great assistance.
>>>>> Many thanks.
>>>>> Ogbos
>>>>>
>>>>>         [[alternative HTML version deleted]]
>>>>>
>>>>> ______________________________________________
>>>>> R-help using r-project.org mailing list -- To UNSUBSCRIBE and more, see
>>>>> https://stat.ethz.ch/mailman/listinfo/r-help
>>>>> PLEASE do read the posting guide http://www.R-project.org/
>>>>> posting-guide.html
>>>>> and provide commented, minimal, self-contained, reproducible code.
>>>>>
>>>>
>>>>
>>

	[[alternative HTML version deleted]]