[R] Monte Carlo on simple regression

[Ogbos Okike]

Dear Erick,

This is great!!
Many thanks for resolving the problem.
Ogbos

On Wed, Aug 22, 2018 at 5:44 PM Eric Berger <ericjberger using gmail.com> wrote:

> 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]]