[R] Monte Carlo simulation for ratio and its CI

Jeff Newmiller jdnewm|| @end|ng |rom dcn@d@v|@@c@@u@
Tue Mar 26 15:26:48 CET 2019


Do you really not know how to use a for loop? The tutorial recommendation seems apropos...

On March 26, 2019 5:57:17 AM PDT, Marna Wagley <marna.wagley using gmail.com> wrote:
>Dear Bert,
>Thank you very much for the response.
>I did it manually but I could not put them in a loop so that I created
>the
>table manually with selecting the rows randomly several times. Here
>what I
>have done so far, please find it. I want to create the table 100 times
>and
>calculate its mean and CI from those 100 values. If anyone can give me
>some
>hint to make a loop, that would be great. I am very grateful with your
>help.
>Thanks,
>
>library(dplyr)
>
>library(plyr)
>
>dat<-structure(list(v1 = c(NA, TRUE, TRUE, TRUE, TRUE, TRUE, NA, TRUE,
>
>NA, NA, TRUE, TRUE, TRUE, TRUE, NA, NA, TRUE, TRUE), v2 = c(TRUE,
>
>NA, NA, NA, NA, TRUE, NA, NA, TRUE, TRUE, NA, TRUE, TRUE, NA,
>
>NA, TRUE, TRUE, NA), v3 = c(TRUE, TRUE, NA, TRUE, TRUE, NA, NA,
>
>TRUE, TRUE, NA, NA, TRUE, TRUE, TRUE, NA, NA, TRUE, NA)), .Names =
>c("v1",
>
>"v2", "v3"), class = "data.frame", row.names = c(NA, -18L))
>
>
>ratio1 <- with(dat, sum(v1,na.rm = TRUE)/sum(v3,na.rm=TRUE))
>
>ratio2 <- with(dat, sum(v2,na.rm = TRUE)/sum(v3,na.rm=TRUE))
>
>#
>
>A1<-sample_n(dat1, 16)# created a table with selecting a 16 sample size
>(rows)
>
>A1.ratio1<-with(A1, sum(v1,na.rm = TRUE)/sum(v3,na.rm=TRUE))
>
>A1.ratio2 <- with(A1, sum(v2,na.rm = TRUE)/sum(v3,na.rm=TRUE))
>
>A1.Table<-data.frame(Ratio1=A1.ratio1, Ratio2=A1.ratio2)
>
>#
>
>A2<-sample_n(dat1, 16)
>
>A2.ratio1<-with(A2, sum(v1,na.rm = TRUE)/sum(v3,na.rm=TRUE))
>
>A2.ratio2 <- with(A2, sum(v2,na.rm = TRUE)/sum(v3,na.rm=TRUE))
>
>A2.Table<-data.frame(Ratio1=A2.ratio1, Ratio2=A2.ratio2)
>
>#
>
>A3<-sample_n(dat1, 16)
>
>A3.ratio1<-with(A3, sum(v1,na.rm = TRUE)/sum(v3,na.rm=TRUE))
>
>A3.ratio2 <- with(A3, sum(v2,na.rm = TRUE)/sum(v3,na.rm=TRUE))
>
>A3.Table<-data.frame(Ratio1=A3.ratio1, Ratio2=A3.ratio2)
>
>#
>
>##..............
>
># I was thinking to repeat this procedure 100 times and calculate the
>ratio
>
>A100<-sample_n(dat1, 16)
>
>A100.ratio1<-with(A100, sum(v1,na.rm = TRUE)/sum(v3,na.rm=TRUE))
>
>A100.ratio2 <- with(A100, sum(v2,na.rm = TRUE)/sum(v3,na.rm=TRUE))
>
>A100.Table<-data.frame(Ratio1=A100.ratio1, Ratio2=A100.ratio2)
>
>#
>
>Tab<-rbind(A1.Table, A2.Table, A3.Table, A100.Table)
>
>
>#Compute the mean for each ratio
>
>Ratio1<-mean(Table1[,1])
>
>Ratio2<-mean(Table1[,2])
>
>
>summary <- ddply(subset(Tab), c(""),summarise,
>
>               N    = length(Tab),
>
>               mean.R1 = mean(Ratio1, na.rm=T),
>
>               median.R1=median(Ratio1, na.rm=T),
>
>               sd.R1   = sd(Ratio1, na.rm=T),
>
>               se.R1   = sd / sqrt(N),
>
>               LCI.95.R1=mean.R1-1.95*se.R1,
>
>               UCI.95.R1=mean.R1+1.95*se.R1,
>
>
>
>               mean.R2 = mean(Ratio2, na.rm=T),
>
>               median.R2=median(Ratio2, na.rm=T),
>
>               sd.R2   = sd(Ratio2, na.rm=T),
>
>               se.R2   = sd / sqrt(N),
>
>               LCI.95.R2=mean.R2-1.95*se.R2,
>
>               UCI.95.R2=mean.R2+1.95*se.R2
>
>               )
>
> summary
>
>
>
>On Mon, Mar 25, 2019 at 4:50 PM Bert Gunter <bgunter.4567 using gmail.com>
>wrote:
>
>>
>> > ratio1 <- with(dat, sum(v1,na.rm = TRUE)/sum(v3,na.rm=TRUE))
>> > ratio1
>> [1] 1.2
>>
>> It looks like you should spend some more time with an R tutorial or
>two.
>> This is basic stuff (if I understand what you wanted correctly).
>>
>> Also, this is not how a "confidence interval" should be calculated,
>but
>> that is another off topic discussion for which
>stats.stackexchange.com is
>> a more appropriate venue.
>>
>> Cheers,
>> Bert
>>
>> Bert Gunter
>>
>> "The trouble with having an open mind is that people keep coming
>along and
>> sticking things into it."
>> -- Opus (aka Berkeley Breathed in his "Bloom County" comic strip )
>>
>>
>> On Mon, Mar 25, 2019 at 4:31 PM Marna Wagley <marna.wagley using gmail.com>
>> wrote:
>>
>>> Hi R User,
>>> I was trying to calculate ratios with confidence interval using
>Monte
>>> Carlo
>>> simulation but I could not figure it out.
>>> Here is the example of my data (see below), I want to calculate
>ratios
>>> (dat$v1/dat$v3 & dat$v2/dat$v3) and its confidence intervals using a
>100
>>> randomly selected data sets.
>>> Could you please give me your suggestions how I can estimate ratios
>with
>>> CI?
>>> I will be very grateful to you.
>>> Sincerely,
>>>
>>> MW
>>> ---
>>> dat<-structure(list(v1 = c(NA, TRUE, TRUE, TRUE, TRUE, TRUE, NA,
>TRUE,
>>>
>>> NA, NA, TRUE, TRUE, TRUE, TRUE, NA, NA, TRUE, TRUE), v2 = c(TRUE,
>>>
>>> NA, NA, NA, NA, TRUE, NA, NA, TRUE, TRUE, NA, TRUE, TRUE, NA,
>>>
>>> NA, TRUE, TRUE, NA), v3 = c(TRUE, TRUE, NA, TRUE, TRUE, NA, NA,
>>>
>>> TRUE, TRUE, NA, NA, TRUE, TRUE, TRUE, NA, NA, TRUE, NA)), .Names =
>c("v1",
>>>
>>> "v2", "v3"), class = "data.frame", row.names = c(NA, -18L))
>>>
>>>
>>> ratio1<-length(which(dat$v1 == "TRUE"))/length(which(dat$v3 ==
>"TRUE"))
>>>
>>> ratio2<-length(which(dat$v2 == "TRUE"))/length(which(dat$v3 ==
>"TRUE"))
>>>
>>>
>>> Thanks
>>>
>>>         [[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]]
>
>______________________________________________
>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.

-- 
Sent from my phone. Please excuse my brevity.



More information about the R-help mailing list