[R] Paternity data analysis problem

John Kane jrkrideau at inbox.com
Wed Jul 24 14:48:51 CEST 2013


Please use dput() to supply data and send in text format not html.

Thanks

John Kane
Kingston ON Canada


> -----Original Message-----
> From: mrahmankufmrt at gmail.com
> Sent: Wed, 24 Jul 2013 19:00:42 +0800
> To: r-help at r-project.org, r-help-request at r-project.org,
> r-help-owner at r-project.org
> Subject: [R] Paternity data analysis problem
> 
> Dear R-helps,
> 
> I did an experiment with FAs ['High' and 'Zero'(no w-3) quality; n=24 for
> each group]. Then I did AI to see their sperm competitiveness based on
> their paternity performance. My data is as below where Fish ID- Blind ID
> for each fish; Group ID- Dietary group ID; Diet quality - High=1, zero=0;
> Babies for paternity- total no. of babies got from females; Success -
> Babies shared/paterned by focal male; Failure - Babies shared/paterned by
> competitor, Proportion - Success/(Success+Failure).
> 
> Fish ID
> 
> Group ID
> 
> Diet quality
> 
> Babies for paternity
> 
> Success
> 
> Failure
> 
> Proportion
> 
> 1
> 
> High
> 
> 1
> 
> 9
> 
> 5
> 
> 4
> 
> 0.556
> 
> 12
> 
> High
> 
> 1
> 
> 7
> 
> 5
> 
> 2
> 
> 0.714
> 
> 15
> 
> High
> 
> 1
> 
> 7
> 
> 4
> 
> 3
> 
> 0.571
> 
> 20
> 
> High
> 
> 1
> 
> 6
> 
> 5
> 
> 1
> 
> 0.833
> 
> 32
> 
> High
> 
> 1
> 
> 7
> 
> 2
> 
> 5
> 
> 0.286
> 
> 37
> 
> High
> 
> 1
> 
> 3
> 
> 1
> 
> 2
> 
> 0.333
> 
> 48
> 
> High
> 
> 1
> 
> 4
> 
> 1
> 
> 3
> 
> 0.25
> 
> 53
> 
> High
> 
> 1
> 
> 10
> 
> 0
> 
> 10
> 
> 0
> 
> 65
> 
> High
> 
> 1
> 
> 3
> 
> 3
> 
> 0
> 
> 1
> 
> 70
> 
> High
> 
> 1
> 
> 4
> 
> 4
> 
> 0
> 
> 1
> 
> 77
> 
> High
> 
> 1
> 
> 7
> 
> 2
> 
> 5
> 
> 0.286
> 
> 82
> 
> High
> 
> 1
> 
> 6
> 
> 6
> 
> 0
> 
> 1
> 
> 96
> 
> High
> 
> 1
> 
> 8
> 
> 2
> 
> 6
> 
> 0.25
> 
> 104
> 
> High
> 
> 1
> 
> 12
> 
> 10
> 
> 2
> 
> 0.833
> 
> 111
> 
> High
> 
> 1
> 
> 4
> 
> 3
> 
> 1
> 
> 0.75
> 
> 123
> 
> High
> 
> 1
> 
> 6
> 
> 5
> 
> 1
> 
> 0.833
> 
> 128
> 
> High
> 
> 1
> 
> 8
> 
> 8
> 
> 0
> 
> 1
> 
> 133
> 
> High
> 
> 1
> 
> 6
> 
> 5
> 
> 1
> 
> 0.833
> 
> 144
> 
> High
> 
> 1
> 
> 12
> 
> 6
> 
> 6
> 
> 0.5
> 
> 152
> 
> High
> 
> 1
> 
> 13
> 
> 11
> 
> 2
> 
> 0.846
> 
> 159
> 
> High
> 
> 1
> 
> 8
> 
> 1
> 
> 7
> 
> 0.125
> 
> 164
> 
> High
> 
> 1
> 
> 4
> 
> 1
> 
> 3
> 
> 0.25
> 
> 169
> 
> High
> 
> 1
> 
> 6
> 
> 2
> 
> 4
> 
> 0.333
> 
> 5
> 
> Zero
> 
> 0
> 
> 9
> 
> 4
> 
> 5
> 
> 0.444
> 
> 10
> 
> Zero
> 
> 0
> 
> 7
> 
> 2
> 
> 5
> 
> 0.286
> 
> 17
> 
> Zero
> 
> 0
> 
> 7
> 
> 3
> 
> 4
> 
> 0.429
> 
> 22
> 
> Zero
> 
> 0
> 
> 6
> 
> 1
> 
> 5
> 
> 0.167
> 
> 36
> 
> Zero
> 
> 0
> 
> 7
> 
> 5
> 
> 2
> 
> 0.714
> 
> 39
> 
> Zero
> 
> 0
> 
> 3
> 
> 2
> 
> 1
> 
> 0.667
> 
> 44
> 
> Zero
> 
> 0
> 
> 4
> 
> 3
> 
> 1
> 
> 0.75
> 
> 51
> 
> Zero
> 
> 0
> 
> 10
> 
> 10
> 
> 0
> 
> 1
> 
> 63
> 
> Zero
> 
> 0
> 
> 3
> 
> 0
> 
> 3
> 
> 0
> 
> 68
> 
> Zero
> 
> 0
> 
> 4
> 
> 0
> 
> 4
> 
> 0
> 
> 73
> 
> Zero
> 
> 0
> 
> 7
> 
> 5
> 
> 2
> 
> 0.714
> 
> 84
> 
> Zero
> 
> 0
> 
> 6
> 
> 0
> 
> 6
> 
> 0
> 
> 94
> 
> Zero
> 
> 0
> 
> 8
> 
> 6
> 
> 2
> 
> 0.75
> 
> 106
> 
> Zero
> 
> 0
> 
> 12
> 
> 2
> 
> 10
> 
> 0.167
> 
> 109
> 
> Zero
> 
> 0
> 
> 4
> 
> 1
> 
> 3
> 
> 0.25
> 
> 121
> 
> Zero
> 
> 0
> 
> 6
> 
> 1
> 
> 5
> 
> 0.167
> 
> 132
> 
> Zero
> 
> 0
> 
> 8
> 
> 0
> 
> 8
> 
> 0
> 
> 137
> 
> Zero
> 
> 0
> 
> 6
> 
> 1
> 
> 5
> 
> 0.167
> 
> 142
> 
> Zero
> 
> 0
> 
> 12
> 
> 6
> 
> 6
> 
> 0.5
> 
> 154
> 
> Zero
> 
> 0
> 
> 13
> 
> 2
> 
> 11
> 
> 0.154
> 
> 157
> 
> Zero
> 
> 0
> 
> 8
> 
> 7
> 
> 1
> 
> 0.875
> 
> 168
> 
> Zero
> 
> 0
> 
> 4
> 
> 3
> 
> 1
> 
> 0.75
> 
> 173
> 
> Zero
> 
> 0
> 
> 6
> 
> 4
> 
> 2
> 
> 0.667
> 
> 
> 
> I ran the following codes to have my results:
> 
> ###Proportion estimate:
> p<-Data$Success/(Data$Success+Data$Failure)
> plot(Data$Group.ID,p,ylab="Proportion of success")
> 
> ###Response variable:
> y<-cbind(Data$Success,Data$Failure)
> model1 <- glm(y~Diet.quality, data=Data, family=binomial)
> summary(model1)
> plot(model1)# gives Q-Q plots
> ###The residual deviance is 152.66  on 44 d.f. so the model is quite
> badly
> overdispersed:
> #152.66/44 where The overdispersion factor is almost 3.46 (unbelievable).
> 
> ## model with logit link functions and weights:
> model2<-glm(cbind(Success,Failure)~Group.ID,data=Data,
> family="binomial"(link="logit"),weights=Success+Failure)
> summary(model2)
> ###The residual deviance is 1196.1  on 46 d.f. so the model is quite
> badly
> overdispersed:
> #1192.1/44 where The overdispersion factor is almost 27.09
> (unbelievable).
> 
> #The simplest way to take this into account is to use what is called an
> #?empirical scale parameter? to reflect the fact that the errors are not
> #binomial as we assumed, but were larger than this (overdispersed) by a
> factor of 3.38.
> 
> model3<-glm(y ~ Group.ID,data=Data,family="quasibinomial")
> summary(model3)
> 
> ###Note that the ratio of the residual deviance and the degrees of
> freedom
> is still
> #larger than 1, but that is no longer a problem as we now allow for
> overdispersion.
> 
>  Each models gives me different results with overdispersion. So, can
> anyone
> help me to give me some valuable suggesions to solve this problem. I'll
> really appreciate your kind assistance and will be grateful to you
> forever.
> 
> With kind regards,
> 
> Moshi
> mrahmankufmrt at gmail.com
> 
> --
> MD. MOSHIUR RAHMAN
> PhD Candidate
> School of Animal Biology/Zoology (M092)
> University of Western Australia
> 35 Stirling Hwy, Crawley, WA, 6009
> Australia.
> Mob.: 061-425205507
> 
> 	[[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.

____________________________________________________________
FREE 3D MARINE AQUARIUM SCREENSAVER - Watch dolphins, sharks & orcas on your desktop!



More information about the R-help mailing list