! Sib pair ACE data
#define nvar 1
#define nvarx2 2
#define nallele 2
#ngroups 2
ACE data
Data NInput=13 ! 13 input variables
missing=-99.989998
Rectangular File=marker6.dat ! read data from file
Labels locus pednum id1 id2 pibd0 pibd1 pibd2 ace1 ace2 a1s1 a2s1 a1s2 a2s2
! aisj is allele i of sibling j
! select only those sib pairs with no missing marker genotypes
select if a1s1 ^= 0
select if a2s1 ^= 0
select if a1s2 ^= 0
select if a2s2 ^= 0
! Phenotype sib 1, phenotype sib 2, probabilities of sharing 0, 1, and 2
! alleles ibd, genotype of sib 1, genotype of sib 2
Select ace1 ace2 pibd0 pibd1 pibd2 a1s1 a2s1 a1s2 a2s2 ;
! use the genotypic information to define the mean and covariance models
Definition_variables pibd0 pibd1 pibd2 a1s1 a2s1 a1s2 a2s2 ;
Begin Matrices;
A Lower nvar nvar Free ! Polygenic component
z lower nvar nvar Free ! Unique environmental component
q lower nvar nvar free ! Linkage component
b full 3 1 ! IBD probabilities (from Genehunter or
! Mapmaker/SIBS)
d full 1 3 ! will contain 0, .5, 1
c fu 1 1 ! constant to contain .5
! association specification
E Full 1 4 ! First allele of sib 1
F Full 1 4 ! Second allele of sib 1
G Full 1 4 ! First allele of sib 2
H Full 1 4 ! Second allele of sib 2
L Full 1 4 ! 0 1 0 1 - matrix to offset to a_w column
m fu 1 1 free ! grand mean
U Unit 4 1 ! Unit matrix
W Full nallele 2 Free ! a_b and a_w parameters, rows are alleles,
! a_b is first column, a_w second column
End Matrices;
! fix the constants
matrix d 0 .5 1
matrix c .5
! starting values for matrices used in linkage component of model
Matrix A .5
Matrix z .7
Matrix q .1
! matrices used for association component of model
! which take genotypes in columns 1 and 3, with columns 2 and 4
! fixed to 1.
Matrix E 1 1 1 1
Matrix F 1 1 1 1
Matrix G 1 1 1 1
Matrix H 1 1 1 1
! L is a constant
Matrix L 0 1 0 1
matrix m 0 ! starting value for estimate of the overall mean
! put ibd probabilities in b
Specify b pibd0 pibd1 pibd2
! put sib genotypes into matrices
Specify E a1s1 0 a1s1 0
Specify F a2s1 0 a2s1 0
Specify G a1s2 0 a1s2 0
Specify H a2s2 0 a2s2 0
Begin Algebra;
! compute the between effects, which are the sum of all the alleles
! in the sibship, obtained by pulling out the appropriate row of
! the first column of the W matrix, which contains the between
! effects for alleles 1 and 2, and summing those effects
I = \part(W,E) + \part(W,F) + \part(W,G) + \part(W,H) ;
! compute the within effects, which use the second column of the W
! matrix. This is the difference between sib 1 and sib 2's
! genotypic scores.
J = \part(W,(E+L)) + \part(W,(F+L)) - \part(W,(G+L)) - \part(W,(H+L)) ;
! compute pihat from ibd probabilities
p = d*b;
End Algebra;
! sib 1's mean is composed of the grand mean plus the between effect
! plus the within effect, sib 2's mean is composed of the grand mean, the
! between effect minus the within effect
Means (M + I + J | M + I - J);
! the expected covariance matrix contains a polygenic component (A),
! a unique environmental component (z) and the QTL component (Q)
Covariance A*A' + z*z' + q*q' | c@a*a' + p@q*q' _
c@a*a' + p@q*q' | a*a' + z*z' + q*q' ;
Option nd=4 ! request 4 decimal places in output
OPtion RS Multiple ! request residuals, multiple fit
Option issat ! this is saturated model for submodel comparison
End
Constrain sum of weights on diffs to equal 0
! this is done because a separate parameter was fitted for each
! of the two possible alleles, while it is only possible to estimate
! one parameter less than the number of alleles
Constraint NI=2
Begin Matrices = Group 1;
O Zero 1 2
P Unit 1 nallele
End Matrices;
! a_b_1 + a_b_2 = 0
! a_w_1 + a_w_2 = 0
Constrain P*W = O ;
Option multiple
option issat ! define the full model
End
! equate a_b to a_w
Specify 1 w
5 5
6 6
option df=-1 ! necessary to get a single df test
End
save combined.mxs
option issat ! a_b = a_w is the new base model
end
! drop both a_b and a_w against a_b = a_w - overall test of association
Drop w 1 1 1 to w 2 2 2
option df=-2
End
! get the a_b = a_w result
! and then drop the linkage parameter
get combined.mxs
drop q 1 1 1 ! test for linkage in the presence of association
end