Hi AK I got the gg plots , thanks. All of a sudden I just got it. I didn't even update packages. Not sure what worked.
Just a basic confusion, to find out if in the association between HiBP and Obese/Overweight status, to quantify if obese males and obese females run a different risk of becoming hypertensive, will the gee output of BP.gee <- gee(HiBP~Obese* Sex ,data=BPsub3,id=CODEA,family=binomial, corstr="exchangeable",na.action=na.omit) suffice? Cheers. > >> On Thu, Jan 17, 2013 at 2:01 PM, arun kirshna [via R] < ml-node+s789695n4655802...@n4.nabble.com> wrote: > HI, > > In the same link at the bottom of the page, > > " > > All is well now after updating all packages with the following: > update.packages()" > > It may or may not solve your problem. > > I got your attachments. You should post those questions in ([hidden > email] <http://user/SendEmail.jtp?type=node&node=4655802&i=0>). I > suggest you to read lme4 book (http://lme4.r-forge.r-project.org/lMMwR/) > #lrgprt.pdf > > A.K. > > > > > > > > ----- Original Message ----- > From: rex2013 <[hidden > email]<http://user/SendEmail.jtp?type=node&node=4655802&i=1>> > > To: [hidden email] <http://user/SendEmail.jtp?type=node&node=4655802&i=2> > Cc: > Sent: Wednesday, January 16, 2013 5:06 AM > Subject: Re: [R] random effects model > > Hi > > I tried removing the missing values and installing "plyr". Still error > message appears with ggplot2 > > Btw, did you get the attachments with my earlier mail? > > Ta. > > On Wed, Jan 16, 2013 at 3:16 AM, arun kirshna [via R] < > [hidden email] <http://user/SendEmail.jtp?type=node&node=4655802&i=3>> > wrote: > > > > > > > Hi, > > Check these links: > > http://comments.gmane.org/gmane.comp.lang.r.ggplot2/6527 > > https://groups.google.com/forum/#!msg/ggplot2/nfVjxL0DXnY/5zf50zCeZuMJ > > A.K. > > > > ________________________________ > > From: Usha Gurunathan <[hidden email]< > http://user/SendEmail.jtp?type=node&node=4655612&i=0>> > > > > To: arun <[hidden email]< > http://user/SendEmail.jtp?type=node&node=4655612&i=1>> > > > > Cc: R help <[hidden email]< > http://user/SendEmail.jtp?type=node&node=4655612&i=2>> > > > > Sent: Tuesday, January 15, 2013 6:31 AM > > Subject: Re: [R] random effects model > > > > > > Hi AK > > > > Got an error message with > > library(ggplot2) > > > > ggplot(BP.stack1,aes(x=factor(HiBP),fill=Obese))+geom_bar(position="fill") > > Error in rename(x, .base_to_ggplot, warn_missing = FALSE) : could not > find > > function "revalue" > > > > ggplot(BP.stack1,aes(x=factor(HiBP),fill=Overweight))+geom_bar(position="fill") > > > Error in rename(x, .base_to_ggplot, warn_missing = FALSE) : could not > find > > function "revalue" > > I got the dot plot, thanks for that. > > > > I have attached some plots, not sure how to interpret, they had unusual > > patterns.Is it because of missing data? I tried removing the missing > data > > too. They still appeared the same. Do I need to transform the data? > > > > > > Thanks in advance. > > > > > > > > > > > > On Tue, Jan 15, 2013 at 8:54 AM, arun <[hidden email]< > http://user/SendEmail.jtp?type=node&node=4655612&i=3>> > > wrote: > > > > HI, > > > > > > > > > > >BP_2b<-read.csv("BP_2b.csv",sep="\t") > > >BP_2bNM<-na.omit(BP_2b) > > > > > >BP.stack3 <- > > > reshape(BP_2bNM,idvar="CODEA",timevar="time",sep="",varying=list(c("Obese14","Obese21"),c("Overweight14","Overweight21"),c("hibp14","hibp21")),v.names=c("Obese","Overweight","HiBP"),times=factor(c(1,2)),direction="long") > > > > > >library(car) > > >BP.stack3$Obese<- recode(BP.stack3$Obese,"1='Obese';0='Not Obese'") > > >BP.stack3$Overweight<- > recode(BP.stack3$Overweight,"1='Overweight';0='Not > > Overweight'") > > > > > >library(ggplot2) > > > >ggplot(BP.stack3,aes(x=factor(HiBP),fill=Obese))+geom_bar(position="fill") > > > > > >ggplot(BP.stack3,aes(x=factor(HiBP),fill=Overweight))+geom_bar(position="fill") > > > > > > > > >You could try lmer() from lme4. > > >library(lme4) > > >fm1<-lmer(HiBP~time+(1|CODEA), family=binomial,data=BP.stack3) #check > > codes, not sure > > >print(dotplot(ranef(fm1,post=TRUE), > > > scales = list(x = list(relation = "free")))[[1]]) > > >qmt1<- qqmath(ranef(fm1, postVar=TRUE)) > > >print(qmt1[[1]]) > > > > > > > > >A.K. > > > > > > > > > > > > > > > > > >________________________________ > > >From: Usha Gurunathan <[hidden email]< > http://user/SendEmail.jtp?type=node&node=4655612&i=4>> > > > > >To: arun <[hidden email]< > http://user/SendEmail.jtp?type=node&node=4655612&i=5>> > > > > >Cc: R help <[hidden email]< > http://user/SendEmail.jtp?type=node&node=4655612&i=6>> > > > > >Sent: Monday, January 14, 2013 6:32 AM > > > > > >Subject: Re: [R] random effects model > > > > > > > > >Hi AK > > > > > >I have been trying to create some plots. All being categorical > variables, > > I am not getting any luck with plots. The few ones that have worked are > > below: > > > > > >barchart(~table(HiBP)|Obese,data=BP.sub3) ## BP.sub3 is the stacked > data > > without missing values > > > > > >barchart(~table(HiBP)|Overweight,data=BP.sub3) > > > > > > >plot(jitter(hibp14,factor=2)~jitter(Obese14,factor=2),col="gray",cex=0.7, > > data=Copy.of.BP_2) ## Copy.of.BP_2 is the original wide format > > > > > >## not producing any good plots with mixed models as well. > > >summary(lme.3 <- lme(HiBP~time, data=BP.sub3,random=~1|CODEA, > > na.action=na.omit)) > > >anova(lme.3) > > >head(ranef(lme.3)) > > >print(plot(ranef(lme.3))) ## > > > > > >Thanks for any help. > > > > > > > > > > > > > > > > > >On Mon, Jan 14, 2013 at 4:33 AM, arun <[hidden email]< > http://user/SendEmail.jtp?type=node&node=4655612&i=7>> > > wrote: > > > > > > > > >> > > >> > > >>HI, > > >> > > >>I think I mentioned to you before that when you reshape the > > >>columns excluding the response variable, response variable gets > repeated > > >>(in this case hibp14 or hibp21) and creates the error" > > >> > > >> > > >>I run your code, there are obvious problems in the code so I didn't > > reach up to BP.gee > > >> > > >> > > >>BP_2b<-read.csv("BP_2b.csv",sep="\t") > > >>BP.stack3 <- > > > reshape(BP_2b,idvar="CODEA",timevar="time",sep="_",varying=list(c("Obese14","Obese21"),c("Overweight14","Overweight21")),v.names=c("Obese","Overweight"),times=factor(c(1,2)),direction="long") > > > > > >> > > >> > > >>BP.stack3 <- > > > transform(BP.stack3,CODEA=factor(CODEA),Sex=factor(Sex,labels=c("Male","Female")),MaternalAge=factor(MaternalAge,labels=c("39years > > > or less","40-49 years","50 years or > > older")),Education=factor(Education,labels=c("Primary/special","Started > > secondary","Completed grade10", "Completed grade12", > > > "College","University")),Birthplace=factor(Birthplace,labels=c("Australia","Other > > > English-speaking","Other"))) > > >> > > >> BP.stack3$Sex <- > > factor(BP.stack3$Sex,levels=levels(BP.stack3$Sex)[c(2,1)]) > > >> > > >> BP.sub3a <- subset(BP.stack3,subset=!(is.na(Sex)| is.na(Education)| > > is.na(Birthplace)|is.na(Education)|is.na(hibp14)| is.na(hibp21))) > > >> nrow(BP.sub3a) > > >>#[1] 3364 > > >> BP.sub5a <- BP.sub3a[order(BP.sub3a$CODEA),] # your code was BP.sub5a > > <- BP.sub3a[order(BP.sub5a$CODEA),] > > >> > > ^^^^^ was not defined before > > >>#Next line > > >>BPsub3$Categ[BPsub6$Overweight==1&BPsub3$time==1&BPsub3$Obese==0]<- > > "Overweight14" #It should be BP.sub3 and what is BPsub6, it was not > > defined previously. > > >>#Error in BPsub3$Categ[BPsub6$Overweight == 1 & BPsub3$time == 1 & > > BPsub3$Obese == : > > >> #object 'BPsub3' not found > > >> > > >> > > >> > > >> > > >> > > >> > > >> > > >>A.K. > > >> > > >> > > >>________________________________ > > >>From: Usha Gurunathan <[hidden email]< > http://user/SendEmail.jtp?type=node&node=4655612&i=8>> > > > > >>To: arun <[hidden email]< > http://user/SendEmail.jtp?type=node&node=4655612&i=9>> > > > > >> > > >>Sent: Sunday, January 13, 2013 1:51 AM > > >> > > >>Subject: Re: [R] random effects model > > >> > > >> > > >> > > >>HI AK > > >> > > >>Thanks a lot for explaining that. > > >> > > >>1. With the chi sq. ( in order to find out if the diffce is > significant > > between groups) do I have create a separate excel file and make a > > dataframe.How do I go about it? > > >> > > >>I have resent a mail to Jun Yan at a difft email ad( first add.didn't > > work, mail not delivered). > > >> > > >>2. With my previous query ( reg. Obese/Overweight/ Normal at age 14 Vs > > change of blood pressure status at 21), even though I had compromised > > without the age-specific regression, but I am still keen to explore why > the > > age-specific regression didn't work despite it looking okay. I have > given > > below the syntax. If you get time, could you kindly look at it and see > if > > it could work by any chance? Apologies for persisting with this query. > > >> > > >> > > >>BP.stack3 <- > > > >>reshape(Copy.of.BP_2,idvar="CODEA",timevar="time",sep="_",varying=list(c("Obese14","Obese21"),c("Overweight14","Overweight21")),v.names=c("Obese","Overweight"),times=factor(c(1,2)),direction="long > > > > > >>BP.stack3 > > >>head(BP.stack3) > > >>tail(BP.stack3) > > >> > > >> names(BP.stack3)[c(2,3,4,5,6,7)] <- > > >>c("Sex","MaternalAge","Education","Birthplace","AggScore","IntScore") > > >> > > >>BP.stack3 <- > > > >>transform(BP.stack3,CODEA=factor(CODEA),Sex=factor(Sex,labels=c("Male","Female")),MaternalAge=factor(MaternalAge,labels=c("39years > > > > > >>or less","40-49 years","50 years or > > > >>older")),Education=factor(Education,labels=c("Primary/special","Started > > >>secondary","Completed grade10", "Completed grade12", > > > >>"College","University")),Birthplace=factor(Birthplace,labels=c("Australia","Other > > > > > >>English-speaking","Other"))) > > >> > > >>table(BP.stack3$Sex) > > >>BP.stack3$Sex <- > > >>factor(BP.stack3$Sex,levels=levels(BP.stack3$Sex)[c(2,1)]) > > >> > > >>levels(BP.stack3$Sex) > > >>BP.sub3a <- subset(BP.stack3,subset=!(is.na(Sex)| is.na(Education)| > > is.na(Birthplace)|is.na(Education)|is.na(hibp14)| is.na(hibp21))) > > >>summary(BP.sub3a) > > >>BP.sub5a <- BP.sub3a[order(BP.sub5a$CODEA),] > > >> BPsub3$Categ[BPsub6$Overweight==1&BPsub3$time==1&BPsub3$Obese==0] > > >><- "Overweight14" > > >>BPsub3$Categ[BPsub6$Overweight==1&BPsub3$time==2&BPsub3$Obese==0] > > >><- "Overweight21" > > > >>BPsub3$Categ[BPsub3$Obese==1&BPsub3$time==1&BPsub3$Overweight==0|BPsub3$Obese==1&BPsub3$time==1&BPsub3$Overweight==1 > > > > > >>] <- "Obese14" > > >>BPsub3$Categ[BPsub3$Obese==0&BPsub3$time==1&BPsub3 > > BPsub3$Categ[BPsub6$Overweight==1&BPsub3$time==1&BPsub3$Obese==0] > > >><- "Overweight14"$Overweight==0] > > >> > > >><- "Normal14" > > >>BPsub3$Categ[BPsub3$Obese==0&BPsub3$time==2&BPsub3$Overweight==0] > > >><- "Normal21" > > > >>BPsub3$Categ[BPsub3$Obese==1&BPsub3$time==2&BPsub3$Overweight==0|BPsub3$Obese==1&BPsub3$time==2&BPsub3$Overweight==1] > > > > > >><- "Obese21" > > >> > > >> > > >> > > >>BPsub3$Categ <- factor(BPsub3$Categ) > > >>BPsub3$time <- factor(BPsub3$time) > > >>summary(BPsub3$Categ) > > >>BPsub7 <- subset(BPsub6,subset=!(is.na(Categ))) > > >>BPsub7$time <- > > >>recode(BPsub7$time,"1=14;2=21") > > >>BPsub7$hibp14 <- factor(BPsub7$hibp14) > > >>levels(BPsub7$hibp14) > > >>levels(BPsub7$Categ) > > >>names(BPsub7) > > >>head(BPsub7) ### this was looking quite okay. > > >> > > >>tail(BPsub7) > > >>str(BPsub7) > > >> > > >>library(gee) > > >> > > >>BP.gee <- geese(hibp14~ time*Categ, > > >>data=BPsub7,id=CODEA,family=binomial, > > >>corstr="exchangeable",na.action=na.omit) > > >> > > >>Thanks again. > > >> > > >> > > >> > > >>On Sun, Jan 13, 2013 at 1:22 PM, arun <[hidden email]< > http://user/SendEmail.jtp?type=node&node=4655612&i=10>> > > wrote: > > >> > > >>HI, > > >>> > > >>>table(BP_2b$Sex) #original dataset > > >>># 1 2 > > >>>#3232 3028 > > >>> nrow(BP_2b) > > >>>#[1] 6898 > > >>> nrow(BP_2bSexNoMV) > > >>>#[1] 6260 > > >>> 6898-6260 > > >>>#[1] 638 #these rows were removed from the BP_2b to create > BP_2bSexNoMV > > >>>BP_2bSexMale<-BP_2bSexNoMV[BP_2bSexNoMV$Sex=="Male",] > > >>> nrow(BP_2bSexMale) > > >>>#[1] 3232 > > >>> nrow(BP_2bSexMale[!complete.cases(BP_2bSexMale),]) #Missing rows > with > > Male > > >>>#[1] 2407 > > >>> nrow(BP_2bSexMale[complete.cases(BP_2bSexMale),]) #Non missing rows > > with Male > > >>>#[1] 825 > > >>> > > >>> > > >>>You did the chisquare test on the new dataset with 6260 rows, right. > > >>>I removed those 638 rows because these doesn't belong to either male > or > > female, but you want the % of missing value per male or female. So, I > > thought this will bias the results. If you want to include the missing > > values, you could do it, but I don't know where you would put that > missing > > values as it cannot be classified as belonging specifically to males or > > females. I hope you understand it. > > >>> > > >>>Sometimes, the maintainer's respond a bit slow. You have to sent an > > email reminding him again. > > >>> > > >>>Regarding the vmv package, you could email Waqas Ahmed Malik ([hidden > > email] <http://user/SendEmail.jtp?type=node&node=4655612&i=11>) > regarding > > options for changing the title and the the font etc. > > >>>You could also use this link ( > > http://www.r-bloggers.com/visualizing-missing-data-2/ ) to plot missing > > value (?plot.missing()). I never used that package, but you could try. > > Looks like it gives more information. > > >>> > > >>>A.K. > > >>> > > >>> > > >>> > > >>> > > >>> > > >>> > > >>> > > >>> > > >>>________________________________ > > >>>From: Usha Gurunathan <[hidden email]< > http://user/SendEmail.jtp?type=node&node=4655612&i=12>> > > > > >>>To: arun <[hidden email]< > http://user/SendEmail.jtp?type=node&node=4655612&i=13>> > > > > >>>Sent: Saturday, January 12, 2013 9:05 PM > > >>> > > >>>Subject: Re: [R] random effects model > > >>> > > >>> > > >>>Hi A.K > > >>> > > >>>So it is number of females missing/total female participants > enrolled: > > 72.65% > > >>>Number of females missing/total (of males+ females) participants > > enrolled : 35.14% > > >>> > > >>>The total no. with the master data: Males: 3232, females: 3028 ( I > got > > this before removing any missing values) > > >>> > > >>>with table(Copy.of.BP_2$ Sex) ## BP > > >>> > > >>> > > >>>If I were to write a table ( and do a chi sq. later), > > >>> > > >>>as Gender Study Non study/missing > > Total > > >>> Male 825 (25.53%) 2407 (74.47%) > > 3232 (100%) > > >>> Female 828 (27.35%) 2200 (72.65%) > 3028 > > ( 100%) > > >>> Total 1653 4607 > > 6260 > > >>> > > >>> > > >>>The problem is when I did > > >>>>colSums(is.na(Copy.of.BP_2), the sex category showed N=638. > > >>> > > >>>I cannot understand the discrepancy.Also, when you have mentioned to > > remove NA, is that not a missing value that needs to be included in the > > total number missing. I am a bit confused. Can you help? > > >>> > > >>>## I tried sending email to gee pack maintainer at the ID with R > site, > > mail didn't go through?? > > >>> > > >>>Many thanks > > >>> > > >>> > > >>> > > >>> > > >>> > > >>> > > >>>On Sun, Jan 13, 2013 at 9:17 AM, arun <[hidden email]< > http://user/SendEmail.jtp?type=node&node=4655612&i=14>> > > wrote: > > >>> > > >>>Hi, > > >>>>Yes, you are right. 72.655222% was those missing among females. > > 35.14377% of values in females are missing from among the whole dataset > > (combined total of Males+Females data after removing the NAs from the > > variable "Sex"). > > >>>> > > >>>>A.K. > > >>>> > > >>>> > > >>>> > > >>>>________________________________ > > >>>>From: Usha Gurunathan <[hidden email]< > http://user/SendEmail.jtp?type=node&node=4655612&i=15>> > > > > >>>>To: arun <[hidden email]< > http://user/SendEmail.jtp?type=node&node=4655612&i=16>> > > > > >>>>Cc: R help <[hidden email]< > http://user/SendEmail.jtp?type=node&node=4655612&i=17>> > > > > >>>>Sent: Saturday, January 12, 2013 5:59 PM > > >>>> > > >>>>Subject: Re: [R] random effects model > > >>>> > > >>>> > > >>>> > > >>>>Hi AK > > >>>>That works. I was trying to get similar results from any other > > package. Being a beginner, I was not sure how to modify the syntax to > get > > my output. > > >>>> > > >>>>lapply(split(BP_2bSexNoMV,BP_ > > >>>>2bSexNoMV$Sex),function(x) > > (nrow(x[!complete.cases(x[,-2]),])/nrow(x))*100) #gives the percentage > of > > rows of missing #values from the overall rows for Males and Females > > >>>>#$Female > > >>>>#[1] 72.65522 > > >>>># > > >>>>#$Male > > >>>>#[1] 74.47401 > > >>>> > > >>>>#iF you want the percentage from the total number rows in Males and > > Females (without NA's in the the Sex column) > > >>>> lapply(split(BP_2bSexNoMV,BP_2bSexNoMV$Sex),function(x) > > (nrow(x[!complete.cases(x[,-2]),])/nrow(BP_2bSexNoMV))*100) > > >>>>#$Female > > >>>>#[1] 35.14377 > > >>>># > > >>>>#$Male > > >>>>#[1] 38.45048 > > >>>> > > >>>>How do I interpret the above 2 difft results? 72.66% of values were > > missing among female participants?? Can you pl. clarify. > > >>>> > > >>>>Many thanks. > > >>>> > > >>>> > > >>>>On Sun, Jan 13, 2013 at 3:28 AM, arun <[hidden email]< > http://user/SendEmail.jtp?type=node&node=4655612&i=18>> > > wrote: > > >>>> > > >>>>lapply(split(BP_2bSexNoMV,BP_2bSexNoMV$Sex),function(x) > > (nrow(x[!complete.cases(x[,-2]),])/nrow(x))*100) #gives the percentage > of > > rows of missing #values from the overall rows for Males and Females > > >>>>>#$Female > > >>>>>#[1] 72.65522 > > >>>>># > > >>>>>#$Male > > >>>>>#[1] 74.47401 > > >>>>> > > >>>>>#iF you want the percentage from the total number rows in Males and > > Females (without NA's in the the Sex column) > > >>>>> lapply(split(BP_2bSexNoMV,BP_2bSexNoMV$Sex),function(x) > > (nrow(x[!complete.cases(x[,-2]),])/nrow(BP_2bSexNoMV))*100) > > >>>>>#$Female > > >>>>>#[1] 35.14377 > > >>>>># > > >>>>>#$Male > > >>>>>#[1] 38.45048 > > >>>> > > >>> > > >> > > > > > > > ______________________________________________ > > [hidden email] > > <http://user/SendEmail.jtp?type=node&node=4655612&i=19>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. > > > > > > ------------------------------ > > If you reply to this email, your message will be added to the > discussion > > below: > > > > > . > > NAML< > http://r.789695.n4.nabble.com/template/NamlServlet.jtp?macro=macro_viewer&id=instant_html%21nabble%3Aemail.naml&base=nabble.naml.namespaces.BasicNamespace-nabble.view.web.template.NabbleNamespace-nabble.view.web.template.NodeNamespace&breadcrumbs=notify_subscribers%21nabble%3Aemail.naml-instant_emails%21nabble%3Aemail.naml-send_instant_email%21nabble%3Aemail.naml> > > > > > > > > -- > View this message in context: > http://r.789695.n4.nabble.com/random-effects-model-tp4654346p4655704.html > > Sent from the R help mailing list archive at Nabble.com. > [[alternative HTML version deleted]] > > ______________________________________________ > [hidden email] <http://user/SendEmail.jtp?type=node&node=4655802&i=4>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. > > > ______________________________________________ > [hidden email] <http://user/SendEmail.jtp?type=node&node=4655802&i=5>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. > > > ------------------------------ > If you reply to this email, your message will be added to the discussion > below: > http://r.789695.n4.nabble.com/random-effects-model-tp4654346p4655802.html > To unsubscribe from random effects model, click > here<http://r.789695.n4.nabble.com/template/NamlServlet.jtp?macro=unsubscribe_by_code&node=4654346&code=dXNoYS5uYXRoYW5AZ21haWwuY29tfDQ2NTQzNDZ8MjAyMjE1NDI0> > . > NAML<http://r.789695.n4.nabble.com/template/NamlServlet.jtp?macro=macro_viewer&id=instant_html%21nabble%3Aemail.naml&base=nabble.naml.namespaces.BasicNamespace-nabble.view.web.template.NabbleNamespace-nabble.view.web.template.NodeNamespace&breadcrumbs=notify_subscribers%21nabble%3Aemail.naml-instant_emails%21nabble%3Aemail.naml-send_instant_email%21nabble%3Aemail.naml> > -- View this message in context: http://r.789695.n4.nabble.com/random-effects-model-tp4654346p4656033.html Sent from the R help mailing list archive at Nabble.com. [[alternative HTML version deleted]] ______________________________________________ R-help@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.