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