! ! QTL via weighted likelihood of three covariance matrices and mean structs ! #define nvar 1 #define startqpar 3 ! nvar*(nvar+1) +1 #define endqpar 3 ! nvar*(nvar+2) #ngroups 6 G1: Calculate variance components and frequency matrices Calculation Begin Matrices; ! Parameters A Full 1 1 Free ! half deviation between homozygotes P Full 1 1 Free ! allele freq D Full 1 1 !Free ! deviation of heterozygote from mean of homozygs E Lower nvar nvar Free ! random envt R Lower nvar nvar Free ! background residual sib resemblance ! Constants H Full 1 1 ! .5 I Unit 1 1 ! for 1-p N Full 1 1 Free ! Grand m mean J Full 3 5 ! IBD 2 coefficients K Full 9 5 ! IBD 1 coefficients L Full 9 5 ! IBD 0 coefficients S Full 1 1 ! 1/16 T Full 1 1 ! 2.0 U unit 3 1 End Matrices; ! Parameters Matrix A 1 ! Matrix D .5 !uncomment for VD Matrix R 1 Matrix E 1 Matrix P .5 Matrix S .0625 Bound .0001 .9999 P 1 1 ! Constants Matrix N 0 Matrix H .5 Matrix T 2 Matrix J ! IBD 2 sixteenths of p^4 p^3q p^2q^2 pq^3 q^4 4 4 4 0 0 ! AA AA 0 4 16 4 0 ! Aa Aa 0 0 4 4 4 ! aa aa Matrix K ! IBD 1 sixteenths of p^4 p^3q p^2q^2 pq^3 q^4 8 4 0 0 0 0 4 8 0 0 0 0 0 0 0 0 4 8 0 0 0 4 16 4 0 0 0 8 4 0 0 0 0 0 0 0 0 8 4 0 0 0 0 4 8 Matrix L ! IBD 0 sixteenths of p^4 p^3q p^2q^2 pq^3 q^4 4 0 0 0 0 0 4 0 0 0 0 0 4 0 0 0 4 0 0 0 0 0 16 0 0 0 0 0 4 0 0 0 4 0 0 0 0 0 4 0 0 0 0 0 4 Begin Algebra; Q = I-P; B = P.P.P.P_P.P.P.Q_P.P.Q.Q_P.Q.Q.Q_Q.Q.Q.Q ; ! p^4 p^3q p^2q^2 pq^3 q^4 V = T.P.Q.(A+D.(Q-P))^T ; ! Additive genetic variance (q^2 in old notation) W = (T.P.Q.D)^T ; ! Dominance genetic variance M = n+A _ n+D _ n-A ; F = M|M ; G = M@U | U@M ; X = J*B@(\sum(j*b)~); ! relative proportions only required so no 1/16 needed Y = K*B@(\sum(k*b)~); Z = L*B@(\sum(l*b)~); End Algebra; End Group; G2:Title build covariance matrices Calculation Begin Matrices; A comp = v1 D comp = w1 E lower 1 1 = E1 R lower 1 1 = R1 H Full 1 1 = H1 U unit 3 1 V unit 9 1 B full 1 1 = P1 J full 1 1 = A1 K full 1 1 = D1 End Matrices; Begin Algebra; S = R.R ; ! Shared residual P = E.E+S; ! Phenotypic variance Q = S; ! Cov of IBD 2 L = S; ! Cov of IBD 1 T = P|Q_Q|P ; ! IBD 2 var-cov O = P|L_L|P ; ! IBD 1 var-cov Z = P|S_S|P ; ! IBD 0 var-cov W = U@T_V@O_V@Z;! Vectorized 21-component covariance matrices X = A|D|E|R|J|K|B; End Algebra; option mxx=x.mat option append option format=(7(F13.6,1x)) End G3:Title Ascertainment Correction for the ASP design Calc Begin Matrices; S Full Free 1 1 !Selection Threshold T Comp = T2 O Comp = O2 Z Comp = Z2 X Unit 1 2 N Full 1 1 =N1 A Full 1 1 =A1 D Full 1 1 =D1 End Matrices; Start 1.28155156554473 S 1 1 ! 10% threshold for N(0,1) Begin Algebra; G = S|S_S|S ; P = N + A; Q = N + D; R = N - A; W = \mnor(T_P|P_G_X)_ \mnor(T_Q|Q_G_X)_ \mnor(T_R|R_G_X)_ \mnor(O_P|P_G_X)_ \mnor(O_P|Q_G_X)_ \mnor(O_P|R_G_X)_ \mnor(O_Q|P_G_X)_ \mnor(O_Q|Q_G_X)_ \mnor(O_Q|R_G_X)_ \mnor(O_R|P_G_X)_ \mnor(O_R|Q_G_X)_ \mnor(O_R|R_G_X)_ \mnor(Z_P|P_G_X)_ \mnor(Z_P|Q_G_X)_ \mnor(Z_P|R_G_X)_ \mnor(Z_Q|P_G_X)_ \mnor(Z_Q|Q_G_X)_ \mnor(Z_Q|R_G_X)_ \mnor(Z_R|P_G_X)_ \mnor(Z_R|Q_G_X)_ \mnor(Z_R|R_G_X) ; End Algebra; End Group; G4:Title QTL model with heterogeneity and weights Data NInput=8 NModel=21 ! READ IN IBD WEIGHTS FOR THE 3 SIBSHIP TYPES RECTANGULAR FILE = qibd.dat LABELS famid marker sib1 sib2 p0 p1 p2 pihat Select if marker = 5 ; Select sib1 sib2 p0 p1 p2 ; Definition p0 p1 p2 ; Begin Matrices ; C Comp = W3 F Comp = F1 G Comp = G1 K Full 1 1 ! p(i)'s go here L Full 1 1 ! p(i)'s go here M Full 1 1 ! p(i)'s go here W Comp = W2 X Comp = X1 Y Comp = Y1 Z Comp = Z1 End Matrices ; Specify K p2 Specify L p1 Specify M p0 Means F_G_G ; Covariance W ; Weights (K@X_L@Y_M@Z)%C ; Options rs samecov multiple issat !func=1e-9 ! option samecov tells it not to re-invert cov matrix as it doesn't change End Group; G5:Title: Mean pi-hat deviation as function of phenotype Data NInput=8 Rectangular file = qibd.dat LABELS famid marker sib1 sib2 p0 p1 p2 pihat Select if marker = 5; Select sib1 sib2 p2 p1 ; Definition sib1 sib2 ; Begin Matrices; S Lower 2 2 free ! variance of pihats X full 1 2 ! observed scores sib1 sib2 T comp = T2 ! IBD 2 var-cov matrix O comp = O2 ! IBD 1 var-cov matrix Z comp = Z2 ! IBD 0 var-cov matrix P comp = P3 ! Mean ibd 2 Q comp = Q3 ! Mean ibd 1 R comp = R3 ! Mean ibd 0 J full 3 5 = J1 K full 9 5 = K1 L full 9 5 = L1 U comp = B1 End Matrices; Matrix S .5 .5 .5 Bound .05 5 S 1 1 to S 2 2 Specify X sib1 sib2 Begin Algebra; A = \pdfnor(X_P|P_T)| ! IBD 2 pdf \pdfnor(X_Q|Q_T)| \pdfnor(X_R|R_T); B = \pdfnor(X_P|P_O)| ! IBD 1 pdf \pdfnor(X_P|Q_O)| \pdfnor(X_P|R_O)| \pdfnor(X_Q|P_O)| \pdfnor(X_Q|Q_O)| \pdfnor(X_Q|R_O)| \pdfnor(X_R|P_O)| \pdfnor(X_R|Q_O)| \pdfnor(X_R|R_O); C = \pdfnor(X_P|P_Z)| ! IBD 0 pdf \pdfnor(X_P|Q_Z)| \pdfnor(X_P|R_Z)| \pdfnor(X_Q|P_Z)| \pdfnor(X_Q|Q_Z)| \pdfnor(X_Q|R_Z)| \pdfnor(X_R|P_Z)| \pdfnor(X_R|Q_Z)| \pdfnor(X_R|R_Z); E = A*(J*U); ! ibd 2 F = B*(K*U); ! ibd 1 G = C*(L*U); ! ibd 0 D = (E)%(E+F+F+G)|(F+F)%(E+F+F+G); ! expected pihat H = S*S'; End Algebra; Mean D; Covariance H; Option RSidual ! Drop @.5 2 End G6: Constrain the ascertainment criteria to 10% of pop. Constraint Begin Matrices; A Comp 1 1 =P3 ! mean of n + a B Comp 1 1 =Q3 ! mean of n + d C Comp 1 1 =R3 ! mean of n - a S Full 1 1 =S3 V Comp = P2 ! variance of trait U Unit 1 1 W Full 1 1 !percent of population selected Q Comp = Q1 D Full 1 1 ! 2 End Matrices; Matrix W .1 Matrix D 2 Begin Algebra; P = (U-Q); K = (P*P)*(\mnor(V_A_S_S_U)) + (D*P*Q)*(\mnor(V_B_S_S_U)) + (Q*Q)*(\mnor(V_C_S_S_U)); End Algebra; Constrain W - K; End ! Eliminate effects of QTL Drop 1 Exit