Skip to content

cumulative sum by group and under some criteria

3 messages · Zjoanna, arun

arun
# 
Hi,
You can also use ?rowMins() or ?rowMaxs() from library(matrixStats)


library(plyr) 
res2<- join(res1,d3,by=c("m1","n1"),type="inner") 
?
p0L<-0.05 
p0H<-0.05 
p1L<-0.20 
p1H<-0.20
?
res2<- within(res2,{p1<- x/m; p2<- y/n;term2_p0<-dbinom(x1,m1, p0L, log=FALSE)* dbinom(y1,n1,p0H, log=FALSE)*dbinom(x-x1,m-m1, p0L, log=FALSE)* dbinom(y-y1,n-n1,p0H, log=FALSE);term2_p1<- dbinom(x1,m1, p1L, log=FALSE)* dbinom(y1,n1,p1H, log=FALSE)*dbinom(x-x1,m-m1, p1L, log=FALSE)* dbinom(y-y1,n-n1,p1H, log=FALSE)})

Pm2<-rbeta(1000, 0.2+res2$x, 0.8+res2$m-res2$x) 
Pn2<-rbeta(1000, 0.2+res2$y, 0.8+res2$n-res2$y) 
Fm2<- ecdf(Pm2) 
Fn2<- ecdf(Pn2) 
library(matrixStats) 
?res3<- within(res2,{Fmm2<-Fm2(p1);Fnn2<- Fn2(p2);R2<- (Fmm2+Fnn2)/2;Fmm_f2<-rowMins(cbind(R2,Fmm2));Fnn_f2<-rowMaxs(cbind(R2,Fnn2));Qm2<- 1-Fmm_f2;Qn2<- 1-Fnn_f2})?
head(res3)
#? m1 n1 x1 y1 m n x y cterm1_P0L cterm1_P1L cterm1_P0H cterm1_P1H?? term2_p1
#1? 2? 2? 0? 0 4 4 0 0???? 0.9025?????? 0.64???? 0.9025?????? 0.64 0.16777216
#2? 2? 2? 0? 0 4 4 0 1???? 0.9025?????? 0.64???? 0.9025?????? 0.64 0.08388608
#3? 2? 2? 0? 0 4 4 0 2???? 0.9025?????? 0.64???? 0.9025?????? 0.64 0.01048576
#4? 2? 2? 0? 0 4 4 1 0???? 0.9025?????? 0.64???? 0.9025?????? 0.64 0.08388608
#5? 2? 2? 0? 0 4 4 1 1???? 0.9025?????? 0.64???? 0.9025?????? 0.64 0.04194304
#6? 2? 2? 0? 0 4 4 1 2???? 0.9025?????? 0.64???? 0.9025?????? 0.64 0.00524288
#????? term2_p0?? p2?? p1?? Qn2 Qm2 Fnn_f2 Fmm_f2???? R2? Fnn2 Fmm2
#1 0.6634204313 0.00 0.00 1.000 1.0? 0.000??? 0.0 0.0000 0.000? 0.0
#2 0.0698337296 0.25 0.00 0.593 1.0? 0.407??? 0.0 0.2035 0.407? 0.0
#3 0.0018377297 0.50 0.00 0.302 1.0? 0.698??? 0.0 0.3490 0.698? 0.0
#4 0.0698337296 0.00 0.25 0.800 0.8? 0.200??? 0.2 0.2000 0.000? 0.4
#5 0.0073509189 0.25 0.25 0.593 0.6? 0.407??? 0.4 0.4035 0.407? 0.4
# 0.0001934452 0.50 0.25 0.302 0.6? 0.698??? 0.4 0.5490 0.698? 0.4
arun
# 
Hello,
?head(d,2)
#? m1 n1 x1 y1 Fmm? Fnn Qm?? Qn? term1_p0 term1_p1
#1? 2? 2? 0? 0?? 0 0.00? 1 1.00 0.8145062?? 0.4096
#2? 2? 2? 0? 1?? 0 0.64? 1 0.36 0.0857375?? 0.2048


Assuming that `x` and `y` are 'x1` and `y1`, where is 'c11` in the dataset?
res2<-do.call(rbind,lapply(lst1[lapply(lst1,nrow)!=0],function(x){
?x[,11:14]<-NA;
?x[,11:12][x$Qm<=c11,]<-cumsum(x[,9:10][x$Qm<=c11,]); 
?x[,13:14][x$Qn<=c12,]<-cumsum(x[,9:10][x$Qn<=c12,]); 
?colnames(x)[11:14]<-
?c("cterm1_P0L","cterm1_P1L","cterm1_P0H","cterm1_P1H");
?x1<-na.locf(x);
?x1[,11:14][is.na(x1[,11:14])]<-0;
?x1}))
Error in `[.data.frame`(x[, 9:10], x$Qm <= c11, ) : 
? object 'c11' not found


A.K.




----- Original Message -----
From: Zjoanna <Zjoanna2013 at gmail.com>
To: r-help at R-project.org
Cc: 
Sent: Sunday, February 24, 2013 4:48 PM
Subject: Re: [R] cumulative sum by group and under some criteria

Thanks! Solved. I have another question.

This is your code for calculated the cumulative sum. how to modify the code
if I want to add another criterion for calculating the cumulative sum:
x[,11:12][x$Qm<=c11,]<-cumsum(x[,9:10][x$Qm<=c11,]);?  # cumsum if
x$Qm<=c11 & x>2
x[,13:14][x$Qn<=c12,]<-cumsum(x[,9:10][x$Qn<=c12,]);? ? #cumsum if x$Qn<=12
& y>2


d <- structure(list(m1 = c(2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,
2, 2, 2, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3),
? ? n1 = c(2, 2, 2, 2, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3, 3, 3, 3,
? ? 3, 3, 3, 3, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2), x1 = c(0,
? ? 0, 0, 1, 1, 1, 2, 2, 2, 0, 0, 0, 0, 1, 1, 1, 1, 2, 2, 2,
? ? 2, 0, 0, 0, 1, 1, 1, 2, 2, 2, 3, 3, 3), y1 = c(0, 1, 2, 0,
? ? 1, 2, 0, 1, 2, 0, 1, 2, 3, 0, 1, 2, 3, 0, 1, 2, 3, 0, 1,
? ? 2, 0, 1, 2, 0, 1, 2, 0, 1, 2), Fmm = c(0, 0, 0, 0.7, 0.59,
? ? 0.64, 1, 1, 1, 0, 0, 0, 0, 0.63, 0.7, 0.74, 0.68, 1, 1, 1,
? ? 1, 0, 0, 0, 0.62, 0.63, 0.6, 0.63, 0.6, 0.68, 1, 1, 1), Fnn = c(0,
? ? 0.64, 1, 0, 0.51, 1, 0, 0.67, 1, 0, 0.62, 0.69, 1, 0, 0.54,
? ? 0.62, 1, 0, 0.63, 0.73, 1, 0, 0.63, 1, 0, 0.7, 1, 0, 0.7,
? ? 1, 0, 0.58, 1), Qm = c(1, 1, 1, 0.65, 0.45, 0.36, 0.5, 0.165,
? ? 0, 1, 1, 1, 1, 0.685, 0.38, 0.32, 0.32, 0.5, 0.185, 0.135,
? ? 0, 1, 1, 1, 0.69, 0.37, 0.4, 0.685, 0.4, 0.32, 0.5, 0.21,
? ? 0), Qn = c(1, 0.36, 0, 0.65, 0.45, 0, 0.5, 0.165, 0, 1, 0.38,
? ? 0.31, 0, 0.685, 0.38, 0.32, 0, 0.5, 0.185, 0.135, 0, 1, 0.37,
? ? 0, 0.69, 0.3, 0, 0.685, 0.3, 0, 0.5, 0.21, 0), term1_p0 =
c(0.81450625,
? ? 0.0857375, 0.00225625, 0.0857375, 0.009025, 0.0002375, 0.00225625,
? ? 0.0002375, 6.25e-06, 0.7737809375, 0.1221759375, 0.00643031249999999,
? ? 0.0001128125, 0.081450625, 0.012860625, 0.000676875, 1.1875e-05,
? ? 0.0021434375, 0.0003384375, 1.78125e-05, 3.125e-07, 0.7737809375,
? ? 0.081450625, 0.0021434375, 0.1221759375, 0.012860625, 0.0003384375,
? ? 0.00643031249999999, 0.000676875, 1.78125e-05, 0.0001128125,
? ? 1.1875e-05, 3.125e-07), term1_p1 = c(0.4096, 0.2048, 0.0256,
? ? 0.2048, 0.1024, 0.0128, 0.0256, 0.0128, 0.0016, 0.32768,
? ? 0.24576, 0.06144, 0.00512, 0.16384, 0.12288, 0.03072, 0.00256,
? ? 0.02048, 0.01536, 0.00384, 0.00032, 0.32768, 0.16384, 0.02048,
? ? 0.24576, 0.12288, 0.01536, 0.06144, 0.03072, 0.00384, 0.00512,
? ? 0.00256, 0.00032)), .Names = c("m1", "n1", "x1", "y1", "Fmm",
"Fnn", "Qm", "Qn", "term1_p0", "term1_p1"), row.names = c(NA,
33L), class = "data.frame")

library(zoo)
lst1<- split(d,list(d$m1,d$n1))
res2<-do.call(rbind,lapply(lst1[lapply(lst1,nrow)!=0],function(x){
x[,11:14]<-NA;
x[,11:12][x$Qm<=c11,]<-cumsum(x[,9:10][x$Qm<=c11,]);?  # cumsum if
x$Qm<=c11 & x>2
x[,13:14][x$Qn<=c12,]<-cumsum(x[,9:10][x$Qn<=c12,]);? ? #cumsum if x$Qn<=12
& y>2
colnames(x)[11:14]<-
c("cterm1_P0L","cterm1_P1L","cterm1_P0H","cterm1_P1H");
x1<-na.locf(x);
x1[,11:14][is.na(x1[,11:14])]<-0;
x1}))
row.names(res2)<- 1:nrow(res2)
res2


On Sat, Feb 23, 2013 at 10:56 PM, arun kirshna [via R] <
ml-node+s789695n4659515h76 at n4.nabble.com> wrote:

            

      
      
      
    

        

      
      
    

        
--
View this message in context: http://r.789695.n4.nabble.com/cumulative-sum-by-group-and-under-some-criteria-tp4657074p4659547.html
Sent from the R help mailing list archive at Nabble.com.
??? [[alternative HTML version deleted]]

______________________________________________
R-help at r-project.org mailing list
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.