m_qtlmap_analyse_lin_gen
[ Top ] [ ANALYSE ] [ Modules ]
NAME
m_qtlmap_analyse_lin_gen
SYNOPSIS
DESCRIPTION
NOTES
SEE ALSO
corniv
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Variables ]
NAME
corniv
DESCRIPTION
get the corresponding estimable index of a level inside a contingence matrix
DIMENSIONS
ntniv
NOTES
contingence
covdir
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Variables ]
NAME
covdir
DESCRIPTION
get the corresponding value of level of a given covariate effect for a kd progeny
DIMENSIONS
nd,ncov
NOTES
le tableau covdir donne pour chaque descendant la valeur des covariables retenues dans le modele contingence
mcov
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Variables ]
NAME
mcov
DESCRIPTION
index of covariate to remove in the construction of the contingence matrix
NOTES
mcov<=0 => no effect, see test_lin contingence
meff
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Variables ]
NAME
meff
DESCRIPTION
index of fixed effect to remove in the construction of the contingence matrix
NOTES
meff<=0 => no effect, see test_lin
mint
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Variables ]
NAME
mint
DESCRIPTION
index of interaction fixed effect-qtl to remove in the construction of the contingence matrix
NOTES
mint<=0 => no effect, see test_lin contingence
nbco
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Variables ]
NAME
nbco
DESCRIPTION
number of covariate in the model for the current trait ic
NOTES
contingence
nbef
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Variables ]
NAME
nbef
DESCRIPTION
number of fixed effect in the model for the current trait ic
NOTES
contingence
nbfem
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Variables ]
NAME
nbfem
DESCRIPTION
number of female estimable
NOTES
contingence
nbniv
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Variables ]
NAME
nbniv
DESCRIPTION
number of level for all fixed effect for the current trait ic
NOTES
contingence
nbnivest
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Variables ]
NAME
nbnivest
DESCRIPTION
number of level estimable (given by the cholesky decomposition and the seuil SEUIL_CHO )
NOTES
nbnivest <= ntniv contingence
nivdir
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Variables ]
NAME
nivdir
DESCRIPTION
get the corresponding level of a given effect for a kd progeny
DIMENSIONS
nd,nteff
NOTES
JM: le tableau nivdir donne pour chaque descendant la position dans le vecteur des effets fixes, des effets exprimes par ce descendant
nteffmax
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Variables ]
NAME
nteffmax
DESCRIPTION
maximum number of effect defined in the contingence matrix
NOTES
contingence
ntlevm
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Variables ]
NAME
ntlevm
DESCRIPTION
number of level for qtl dam interaction with fixed effect
NOTES
contingence
ntlevp
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Variables ]
NAME
ntlevp
DESCRIPTION
number of level for qtl sire interaction with fixed effect
NOTES
contingence
ntnifix
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Variables ]
NAME
ntnifix
DESCRIPTION
index of the latest level of fixed effect (use to get first level/effect of covariate)
NOTES
contingence
ntniv
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Variables ]
NAME
ntniv
DESCRIPTION
number of level defined in the contingence matrix
NOTES
contingence
ntnivmax
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Variables ]
NAME
ntnivmax
DESCRIPTION
maximum number of level (maximum column size) in the contingence matrix
NOTES
contingence
par0
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Variables ]
NAME
par0
DESCRIPTION
Solution finded under H0
DIMENSIONS
nbnivest
NOTES
par1
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Variables ]
NAME
par1
DESCRIPTION
Solution finded under H1
DIMENSIONS
nbnivest
NOTES
pb_haplo
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Variables ]
NAME
pb_haplo
DESCRIPTION
DIMENSIONS
nd,NB_HAPLO_PRIOR
NOTES
preprinc,contingence_ldla
pm_ldla
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Variables ]
NAME
pm_ldla
DESCRIPTION
probabilities to receive the ith haplotype dam for the kkd progenies according the genotype, kkd come from ngend
DIMENSIONS
maxval(ngenom),maxval(ngend),2
NOTES
preprinc,contingence_ldla
pmt
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Variables ]
NAME
pmt
DESCRIPTION
DIMENSIONS
NOTES
pp_ldla
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Variables ]
NAME
pp_ldla
DESCRIPTION
probabilities to receive the ith haplotype sire for the kkd progenies, kkd come from ngend
DIMENSIONS
maxval(ngend),2
NOTES
preprinc,contingence_ldla
ppt
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Variables ]
NAME
ppt
DESCRIPTION
DIMENSIONS
NOTES
prbm
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Variables ]
NAME
prbm
DESCRIPTION
probabilities cumulates (according all dam genotype probabilities) to receive the qtl from dam
DIMENSIONS
nd
NOTES
preprinc,contingence
prbp
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Variables ]
NAME
prbp
DESCRIPTION
probabilities cumulates (according all dam genotype probabilities) to receive the qtl from sire
DIMENSIONS
nd
NOTES
preprinc,contingence
precis
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Variables ]
NAME
precis
DESCRIPTION
precision of each level
DIMENSIONS
ntniv
NOTES
contingence
precis0
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Variables ]
NAME
precis0
DESCRIPTION
Precision of each level under H0
DIMENSIONS
ntniv
NOTES
precis1
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Variables ]
NAME
precis1
DESCRIPTION
Precision of each level under H1
DIMENSIONS
ntniv
NOTES
prop_haplo_info
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Variables ]
NAME
prop_haplo_info
DESCRIPTION
DIMENSIONS
NOTES
vecsol
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Variables ]
NAME
vecsol
DESCRIPTION
boolean vector of estimability of level referenced in the contingence matrix
DIMENSIONS
ntniv
NOTES
contingence
vecsol0
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Variables ]
NAME
vecsol0
DESCRIPTION
boolean vector of estimability level under H0
DIMENSIONS
ntniv
NOTES
vecsol1
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Variables ]
NAME
vecsol1
DESCRIPTION
boolean vector of estimability level under H1
DIMENSIONS
ntniv
NOTES
xinc
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Variables ]
NAME
xinc
DESCRIPTION
contingence matrix
DIMENSIONS
nd,ntniv
NOTES
xx
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Variables ]
NAME
xx
DESCRIPTION
the incidence matrix X'.X , X contingence matrix
DIMENSIONS
ntniv,ntniv
NOTES
contingence
xxx
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Variables ]
NAME
xxx
DESCRIPTION
the inverse of the incidence matrix (X'.X) -1 , X contingence matrix
DIMENSIONS
ntniv,ntniv
NOTES
contingence,confusion
confusion
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Subroutines ]
NAME
confusion
DESCRIPTION
Print in the result file, the confusion (real values) between qtl effet and the other effect
INPUTS
mod : 'avant','apres'
ic : index of the trait
chr : index of the chromosome
est_moy : true => the general mean is included as an effect
option_anal : type analysis * 'LA ', 'LD ', 'LDLA', 'LDJH' *
NOTES
Mesure de la confuqion entre les effets QTL et d'autres effets
THRES_CONFUSION
SOURCE
1814 subroutine confusion(mod,chr,ic,est_moy,option_anal) 1815 1816 integer , intent(in) :: ic,chr 1817 character(len=*) , intent(in) :: mod 1818 logical , intent(in),optional :: est_moy 1819 character(len=4) , intent(in) :: option_anal 1820 real (kind=dp) :: xcor(ntnivmax,ntnivmax) 1821 logical , dimension(nfem) :: femIsOk 1822 ! real (kind=dp) ,dimension(ntniv,ntniv) :: xxx 1823 ! 1824 ! Divers 1825 logical :: icas,itest 1826 integer :: nfirst,ief,i,j,jef,ip,ilevp,iclic,iqtl,ntot 1827 integer :: jniv,ilev,iniveau,ico,jm,ilevm,ifem 1828 real (kind=dp) ::corrmin 1829 logical :: estMoyLocal = .true. 1830 ! 1831 1832 if ( present(est_moy) ) then 1833 estMoyLocal = est_moy 1834 end if 1835 ! 1836 if(mod.eq.'avant') then 1837 itest=.false. 1838 call prepinc(chr,1,ic,option_anal) 1839 call contingence(ic,1,itest,estMoyLocal) 1840 write(nficout,601) 1841 601 format(//,80('*')/'Test of confusion between QTL ', & 1842 'and other effects in the initial full model',/, & 1843 '(test based on the correlation between columns of the ','incidence matrix)',//) 1844 end if 1845 1846 if(mod.eq.'apres') then 1847 write(nficout,602) 1848 602 format(//,80('*')/'Test of confusion between QTL ', & 1849 'and other effects in the final constained model',/, & 1850 '(test based on the correlation between columns of the ','incidence matrix)',//) 1851 end if 1852 ! icas est un indicateur de defaut possible d'identifiabilit� 1853 ! 1854 icas=.false. 1855 ! 1856 1857 if (estMoyLocal) then 1858 nfirst= 1+ntlevp*np+ntlevm*namest(ic)+1 1859 else 1860 nfirst= ntlevp*np+ntlevm*namest(ic)+1 1861 end if 1862 1863 corrmin = 0.d0 1864 1865 write(nficout,600) 1866 600 format(//,80('*')/'Confusion between QTL and other effects ','(final constained model)',//) 1867 1868 ! 1869 ! xcor(i,j), la correlation entre les estim�es i et j d'apr�s X'X-1 1870 ! 1871 xcor=0.d0 1872 ief=0 1873 do i=1,ntniv 1874 if(vecsol(i))then 1875 ief=ief+1 1876 jef=0 1877 do j=1,ntniv 1878 if(vecsol(j))then 1879 jef=jef+1 1880 if (xxx(ief,ief) == 0 .or. xxx(jef,jef) == 0) then 1881 xcor(i,j)=0.d0 1882 else 1883 xcor(i,j)=xxx(ief,jef)/dsqrt(xxx(ief,ief)*xxx(jef,jef)) 1884 end if 1885 end if 1886 end do 1887 end if 1888 end do 1889 1890 ! 1891 ! on caclule les correlations entre les colonnes de X'X des effets qtl p�re et de 1892 ! tous les autres effets 1893 ! 1894 do ip=1,np 1895 ! print *,"ip:",ip 1896 do ilevp=1,ntlevp 1897 iclic=0 1898 if (estMoyLocal) then 1899 iqtl=1+ilevp+ntlevp*(ip-1) 1900 else 1901 iqtl=ilevp+ntlevp*(ip-1) 1902 end if 1903 1904 if (vecsol(iqtl))then 1905 ! 1906 ! test d'une confusion avec les effets polyg�niques p�re 1907 ! 1908 ntot=nfirst-1 1909 !print *,"iqtl:",iqtl 1910 do ief=1,np 1911 jniv=ntot+ief 1912 ! print *,'jniv:',jniv 1913 ! print *,'jniv:',jniv,dabs(xcor(iqtl,jniv)) 1914 if(dabs(xcor(iqtl,jniv)).ge.corrmin)corrmin=dabs(xcor(iqtl,jniv)) 1915 ! 1916 ! si la correlation est sup�rieure � corMax, on imprime une alerte 1917 ! 1918 if(dabs(xcor(iqtl,jniv)).ge.THRES_CONFUSION)then 1919 icas=.true. 1920 if(iclic.eq.0) then 1921 write(nficout,605)trim(pere(ip)),ilevp 1922 605 format('Risk for sire ',associated,' of confusion between the QTL level',i3,' and :') 1923 iclic=1 1924 end if 1925 write(nficout,610)trim(pere(ief)),xcor(iqtl,jniv) 1926 610 format('the sire ',associated,' polygenic effect ',' (correlation ',f5.2,')') 1927 end if 1928 end do 1929 1930 ! 1931 ! test d'une confusion avec les effets polyg�niques m�re 1932 ! 1933 ntot=ntot+np 1934 1935 do ief=1,nbfem 1936 jniv=ntot+ief 1937 ! 1938 ! si la correlation est sup�rieur � seuil, on imprime une alerte 1939 ! 1940 if(dabs(xcor(iqtl,jniv)).ge.corrmin)corrmin=dabs(xcor(iqtl,jniv)) 1941 if(dabs(xcor(iqtl,jniv)).ge.THRES_CONFUSION)then 1942 icas=.true. 1943 if(iclic.eq.0) then 1944 write(nficout,605)trim(pere(ip)),ilevp 1945 iclic=1 1946 end if 1947 write(nficout,611)trim(mere(ief)),xcor(iqtl,jniv) 1948 611 format('the dam ',associated,' polygenic effect ',' (correlation ',f5.2,')') 1949 end if 1950 end do 1951 ! 1952 ! test d'une confusion avec les effets fix�s 1953 ! 1954 1955 ntot=ntot+nbfem 1956 if(nbniv.ne.0)then 1957 ilev=0 1958 do ief=1,nbef 1959 do iniveau=1,nlev(modele(ic,3+ief)) 1960 ilev=ilev+1 1961 jniv=ntot+ilev 1962 ! 1963 ! si la correlation est sup�rieur � corMax, on imprime une alerte 1964 ! 1965 if(dabs(xcor(iqtl,jniv)).ge.corrmin)corrmin=dabs(xcor(iqtl,jniv)) 1966 if(dabs(xcor(iqtl,jniv)).ge.THRES_CONFUSION)then 1967 icas=.true. 1968 if(iclic.eq.0) then 1969 write(nficout,605)trim(pere(ip)),ilevp 1970 iclic=1 1971 end if 1972 1973 write(nficout,612)iniveau,trim(namefix(modele(ic,3+ief))),xcor(iqtl,jniv) 1974 612 format('the level ',i3,' of the effect ',a15,' (correlation ',f5.2,')') 1975 end if 1976 end do 1977 end do 1978 end if 1979 ! 1980 ! test d'une confusion avec les covariables 1981 ! 1982 ntot=ntot+nbniv 1983 if(nbco.ne.0)then 1984 do ico=1,nbco 1985 jniv=ntot+ico 1986 ! 1987 ! si la correlation est sup�rieur � corMax, on imprime une alerte 1988 ! 1989 if(dabs(xcor(iqtl,jniv)).ge.corrmin)corrmin=dabs(xcor(iqtl,jniv)) 1990 if(dabs(xcor(iqtl,jniv)).ge.THRES_CONFUSION)then 1991 icas=.true. 1992 if(iclic.eq.0) then 1993 write(nficout,605)trim(pere(ip)),ilevp 1994 iclic=1 1995 end if 1996 write(nficout,613)trim(namecov(modele(ic,3+nbef+ico))),xcor(iqtl,jniv) 1997 613 format('the covariable ',a15,' (correlation ',f5.2,')') 1998 end if 1999 end do 2000 ! 2001 end if 2002 end if 2003 ! 2004 ! fin pour le niveau du qtl du pere 2005 ! 2006 end do 2007 end do 2008 2009 ! 2010 ! meme scenario pour les m�res 2011 ! 2012 femisok=.false. 2013 do jm=1,nm 2014 if (.not. estime(ic,jm)) cycle 2015 ifem=iam(ic,repfem(jm)) 2016 if (femIsok(ifem)) cycle 2017 femisok(ifem)=.true. 2018 ! print *,"jm:",jm 2019 do ilevm=1,ntlevm 2020 iclic=0 2021 2022 if (estMoyLocal) then 2023 iqtl=1+ntlevp*np+ilevm+ntlevm*(ifem-1) 2024 else 2025 iqtl=ntlevp*np+ilevm+ntlevm*(ifem-1) 2026 end if 2027 ! print *,"iqtl:",iqtl 2028 if (vecsol(iqtl))then 2029 ! 2030 ! test d'une confusion avec les effets polyg�niques p�re 2031 ! 2032 ntot=nfirst-1 2033 do ief=1,np 2034 jniv=ntot+ief 2035 ! 2036 ! si la correlation est sup�rieure � corMax, on imprime une alerte 2037 ! 2038 ! print *,'jniv:',jniv 2039 ! print *,dabs(xcor(iqtl,jniv)) 2040 if(dabs(xcor(iqtl,jniv)).ge.THRES_CONFUSION)then 2041 icas=.true. 2042 if(iclic.eq.0) then 2043 write(nficout,615)trim(mere(jm)),ilevm 2044 615 format('Risk for dam ',associated,' of confusion between the QTL level',i3,' and :') 2045 iclic=1 2046 end if 2047 write(nficout,610)trim(pere(ief)),xcor(iqtl,jniv) 2048 end if 2049 end do 2050 ! 2051 ! test d'une confusion avec les effets polyg�niques m�re 2052 ! 2053 ntot=ntot+np 2054 do ief=1,nbfem 2055 jniv=ntot+ief 2056 ! 2057 ! si la correlation est sup�rieur � corMax, on imprime une alerte 2058 ! 2059 if(dabs(xcor(iqtl,jniv)).ge.corrmin)corrmin=dabs(xcor(iqtl,jniv)) 2060 if(dabs(xcor(iqtl,jniv)).ge.THRES_CONFUSION)then 2061 icas=.true. 2062 if(iclic.eq.0) then 2063 write(nficout,615)trim(mere(jm)),ilevm 2064 iclic=1 2065 end if 2066 write(nficout,611)trim(mere(ief)),xcor(iqtl,jniv) 2067 end if 2068 end do 2069 ! 2070 ! test d'une confusion avec les effets fix�s 2071 ! 2072 2073 ntot=ntot+nbfem 2074 if(nbniv.ne.0)then 2075 ilev=0 2076 do ief=1,nbef 2077 do iniveau=1,nlev(modele(ic,3+ief)) 2078 ilev=ilev+1 2079 jniv=ntot+ilev 2080 ! 2081 ! si la correlation est sup�rieur � corMax, on imprime une alerte 2082 ! 2083 if(dabs(xcor(iqtl,jniv)).ge.corrmin)corrmin=dabs(xcor(iqtl,jniv)) 2084 if(dabs(xcor(iqtl,jniv)).ge.THRES_CONFUSION)then 2085 icas=.true. 2086 if(iclic.eq.0) then 2087 write(nficout,615)trim(mere(jm)),ilevm 2088 iclic=1 2089 end if 2090 2091 write(nficout,612)iniveau,trim(namefix(modele(ic,3+ief))),xcor(iqtl,jniv) 2092 end if 2093 end do 2094 end do 2095 end if 2096 ! 2097 ! test d'une confusion avec les covariables 2098 ! 2099 ntot=ntot+nbniv 2100 if(nbco.ne.0)then 2101 do ico=1,nbco 2102 jniv=ntot+ico 2103 ! 2104 ! si la correlation est sup�rieur � corMax, on imprime une alerte 2105 ! 2106 if(dabs(xcor(iqtl,jniv)).ge.corrmin)corrmin=dabs(xcor(iqtl,jniv)) 2107 if(dabs(xcor(iqtl,jniv)).ge.THRES_CONFUSION)then 2108 icas=.true. 2109 if(iclic.eq.0) then 2110 write(nficout,615)trim(mere(jm)),ilevm 2111 iclic=1 2112 end if 2113 write(nficout,613)trim(namecov(modele(ic,3+nbef+ico))),xcor(iqtl,jniv) 2114 end if 2115 end do 2116 ! 2117 end if 2118 end if 2119 ! 2120 2121 ! fin pour le niveau du qtl de la mere 2122 ! 2123 end do 2124 end do 2125 ! 2126 ! situation sans probleme 2127 ! 2128 if(.not.icas) then 2129 write(nficout,616)corrmin 2130 616 format(/,' No confusion detected',/, & 2131 ' the highest correlation is : ', F7.3,/,80('*')/) 2132 end if 2133 2134 return 2135 end subroutine confusion
contingence
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Subroutines ]
NAME
contingence
DESCRIPTION
INPUTS
ic : index of the trait
ihyp : under this hypothesis
itest : testing nuisances effects
est_moy : estimation of the general mean
NOTES
construction de la matrice de contingence et recherche des contraintes
SOURCE
625 subroutine contingence (ic,ihyp,itest,est_moy) 626 ! 627 ! Tableaux dimensionn�s au nombre d'effets � estimer 628 integer :: indestim(nm) 629 ! 630 ! divers 631 integer i,j,k,iniv,jniv,ki,kj,kniv,kkd,kd,ief,jef,icov,jcov 632 integer km,ip,nm1,nm2,nd1,nd2,jm,nbtef,nbtco,nbint,nbtint 633 integer nkd,nteff,nivintp,nbt,nivintm,nd 634 logical , intent(in) :: itest 635 integer , intent(in) :: ic,ihyp 636 logical , intent(in) , optional :: est_moy 637 real (kind=dp) , dimension(:,:),allocatable :: triang 638 integer :: stat 639 640 logical :: estmoylocal=.true. 641 642 ! NOTE Olivier : j ai mis est_moy optional et j initlise une variable local 643 ! modlin de jm n a pas besoin d etre modifie, c est donc moin intrusif.... 644 if ( present(est_moy) ) then 645 estmoylocal=est_moy 646 end if 647 648 !cccccccccccccccccccccccc 649 !cccccccccccccccccccccccc 650 ! L'analyse est faite pour le ic �me caract�re 651 ! 652 ! Les premier �l�ment du vecteur des effets � estimer sont la moyenne 653 ! g�n�rale, les effets p�re et les effets m�re 654 ! 655 ! 656 ! initilisation des compteurs 657 ! 658 ! nbfem est le nombre de m�res "estimables" d'apr�s ndmin 659 ! nbef est le nombre d'effets fix�s parasites 660 ! nbco le nombre de covariables 661 ! nbintp le nombre d'effet hierarchis�s dnas les effets qtl p�re 662 ! nbintm idem pour les effets qtl m�re 663 ! nbniv est le nombre total de niveaux d'effets parasites 664 ! nbnivp, celui des efets hierarchis�s dans les effets qtl p�re 665 ! nbnivm, idem pour les qtl m�re 666 ! 667 ! 668 nbfem=0 669 do jm=1,nm 670 if(estime(ic,jm)) then 671 nbfem=nbfem+1 672 indestim(jm)=nbfem 673 end if 674 end do 675 nbef=modele(ic,1) 676 nbtef=nbef 677 nbco=modele(ic,2) 678 nbtco=nbco 679 nbint=modele(ic,3) 680 nbtint=nbint 681 ! 682 ! si itest est � true, il s'agit de retirer un des effet du modele 683 ! et de recalculer la matrice X'X et les estimabilit� 684 ! 685 if(.not.itest) then 686 meff=0 687 mcov=0 688 mint=0 689 else 690 if(meff.ne.0) nbef=nbef-1 691 if(mcov.ne.0) nbco=nbco-1 692 if(mint.ne.0) nbint=nbint-1 693 end if 694 695 ! 696 ! construction de la matrice d' incidence xinc 697 ! de dimension nb de descendants x 1+nb pere + nb mere + 698 ! nb niveaux effets fix�s + nb c covariables (sous H0) 699 ! 700 ! 701 ! le tableau nivdir donne pour chaque descendant la position 702 ! dans le vecteur des effets fix�s, des effets exprim�s par ce descendant 703 ! 704 ! le tableau covdir donne pour chaque descendant la valeur des covariables 705 ! retenues dans le mod�le 706 ! 707 ! ntniv est le nombre total de niveaux d'effets et covariables 708 ! 709 nkd=0 710 xinc=0.d0 711 nivdir = 0 712 do ip=1,np 713 nm1=nmp(ip)+1 714 nm2=nmp(ip+1) 715 do jm=nm1,nm2 716 nd1=ndm(jm)+1 717 nd2=ndm(jm+1) 718 do kd=nd1,nd2 719 if(presentc(ic,kd)) then 720 nkd=nkd+1 721 ! 722 ! le premier effet fix� est la moyenne g�n�rale 723 ! 724 if (estmoylocal) then 725 nivdir(kd,1)=1 726 xinc(nkd,1)=1 727 ntniv=1 728 nteff=1 729 else ! pas d esitmation de la moyenne generale 730 ntniv=0 731 nteff=0 732 endif 733 734 ! 735 ! puis, sous H1, les effets QTL 736 ! 737 if(ihyp.eq.1) then 738 ! 739 ! si des effets fix�s sont en interaction avec l'all�le paternel au qtl 740 ! il faut estimer les effets qtl pour chacun des niveaux concern�s 741 ! 742 ! on commence donc par cr�er un effet composites regroupant l'ensemble de 743 ! ces effets 744 ! 745 nivintp=1 746 ntlevp=1 747 ! 07/04/2010 748 ! correction bug ofi : nbt=3+nbef+nbco n est pas valide 749 nbt=3+nbtef+nbtco 750 if(nbint.ne.0) then 751 ief=0 752 do jef=1,nbtint 753 if(mint.ne.jef)then 754 ief=ief+1 755 nivintp=nivintp+ntlevp*(nivx(kd,modele(ic,nbt+jef))-1) 756 ntlevp=ntlevp*nlev(modele(ic,nbt+jef)) 757 end if 758 end do 759 end if 760 761 ! 762 ! puis on cr�e le tableau des niveaux cumul�s et la matrice d'incidence 763 ! 764 nivdir(kd,nteff+1)=ntniv+nivintp+ntlevp*(ip-1) 765 xinc(nkd,nivdir(kd,nteff+1))=prbp(kd) 766 ntniv=ntniv+ntlevp*np 767 nteff=nteff+1 768 ! 769 ! m�me op�ration pour le qtl transmis par la m�re 770 ! 771 if(nbfem.ne.0)then 772 nivintm=1 773 ntlevm=1 774 ! 07/04/2010 775 ! correction bug ofi : nbt=nbt+nbint n est pas valide 776 !nbt=nbt+nbint 777 if(nbint.ne.0) then 778 ief=0 779 do jef=1,nbtint 780 if(nbint.ne.0) then 781 ief=ief+1 782 nivintm=nivintm+ntlevm*(nivx(kd,modele(ic,nbt+jef))-1) 783 ntlevm=ntlevm*nlev(modele(ic,nbt+jef)) 784 end if 785 end do 786 end if 787 788 if(estime(ic,jm)) then 789 km=iam(ic,repfem(jm)) 790 nivdir(kd,nteff+1)=ntniv+nivintm+ntlevm*(km-1) 791 xinc(nkd,nivdir(kd,nteff+1))=prbm(kd) 792 end if 793 ntniv=ntniv+ntlevm*namest(ic) 794 nteff=nteff+1 795 end if 796 ! 797 ! fin de H1 798 ! 799 end if 800 ! 801 ! puis les effets polyg�niques parentaux 802 ! 803 nivdir(kd,nteff+1)= ntniv+ip 804 xinc(nkd,nivdir(kd,nteff+1))=1 805 ntniv=ntniv+np 806 nteff=nteff+1 807 808 if(nbfem.ne.0)then 809 if(estime(ic,jm)) then 810 nivdir(kd,nteff+1)=ntniv+indestim(jm) 811 xinc(nkd,nivdir(kd,nteff+1))=1 812 end if 813 ntniv=ntniv+nbfem 814 nteff=nteff+1 815 end if 816 817 ! 818 ! autres effets fix�s 819 ! 820 nbniv=0 821 ief=0 822 do jef=1,nbtef 823 if(meff.ne.jef) then 824 ief=ief+1 825 nivdir(kd,nteff+ief)=ntniv+nbniv+nivx(kd,modele(ic,3+jef)) 826 xinc(nkd,nivdir(kd,nteff+ief))=1 827 nbniv=nbniv+nlev(modele(ic,3+jef)) 828 end if 829 end do 830 831 ntnifix=ntniv+nbniv 832 ! 07/04/2010 833 ! correction bug ofi : nbt=3+nbef n est pas valide 834 nbt=3+nbtef 835 836 ! 837 ! covariables 838 ! 839 icov=0 840 do jcov=1,nbtco 841 if(mcov.ne.jcov) then 842 icov=icov+1 843 covdir(kd,icov)=covar(kd,modele(ic,nbt+jcov)) 844 xinc(nkd,ntnifix+icov)=covar(kd,modele(ic,nbt+jcov)) 845 end if 846 end do 847 ntniv=ntnifix+nbco 848 nbt=nbt+nbco 849 850 end if 851 852 853 end do 854 end do 855 end do 856 857 allocate (triang(ntniv,ntniv),STAT=stat) 858 call check_allocate(stat,'triang [m_qtlmap_analyse_lin_gen]') 859 860 ! 861 ! construction de X'X 862 ! 863 xx=0.d0 864 865 do kd=1,nkd 866 do iniv =1,ntniv 867 do jniv=1,ntniv 868 xx(iniv,jniv)=xx(iniv,jniv)+xinc(kd,iniv)*xinc(kd,jniv) 869 end do 870 end do 871 end do 872 873 ! application de la d�composition de choleski pour d�terminer les contraintes 874 ! 875 ! Dans le vecteur vecsol les effets � estimer sont indiqu�s � TRUE et ceux 876 ! qui ne sont pas estimables � FALSE 877 ! 878 ! La matrice triang est la matrice triangulaire sup�rieure M telle que 879 ! M'M = X'X 880 ! 881 do i=1,ntniv 882 do j=1,ntniv 883 triang(i,j)=0.d0 884 end do 885 end do 886 887 do j=1,ntniv 888 vecsol(j)=.true. 889 triang(j,j)=xx(j,j) 890 do k=1,j-1 891 triang(j,j)=triang(j,j)- triang(j,k)*triang(j,k) 892 end do 893 894 if(triang(j,j).gt.SEUIL_CHO) then 895 triang(j,j)=dsqrt(triang(j,j)) 896 do i=j+1,ntniv 897 triang(i,j)=xx(i,j) 898 do k=1,j-1 899 triang(i,j)=triang(i,j)-triang(i,k)*triang(j,k) 900 end do 901 triang(i,j)=triang(i,j)/triang(j,j) 902 end do 903 904 else 905 vecsol(j)=.false. 906 end if 907 end do 908 909 ! 910 ! Table de correspondance entre niveaux de d�part et niveaux estimables 911 ! Au lieu de ntniv effets � estimer, on en a plus que nbnivest 912 ! corniv(i) est la position, dans la liste des effets estimables 913 ! du i �me effet initial 914 ! 915 916 nbnivest=0 917 do ief=1,ntniv 918 if(vecsol(ief))then 919 nbnivest=nbnivest+1 920 corniv(ief)=nbnivest 921 end if 922 end do 923 924 return 925 end subroutine contingence
contingence_cox
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Subroutines ]
NAME
contingence_cox
DESCRIPTION
INPUTS
ic : index of the trait
ihyp : under this hypothesis
itest : testing nuisances effects
est_moy : estimation of the general mean
NOTES
construction de la matrice de contingence et recherche des contraintes
SOURCE
942 subroutine contingence_cox (ic,ihyp,itest,est_moy) 943 ! 944 ! Tableaux dimensionn�s au nombre d'effets � estimer 945 integer :: indestim(nm) 946 ! 947 ! divers 948 integer i,j,k,iniv,jniv,ki,kj,kniv,kkd,kd,ief,jef,icov,jcov 949 integer km,ip,nm1,nm2,nd1,nd2,jm,nbtef,nbtco,nbint,nbtint 950 integer nkd,nteff,nivintp,nbt,nivintm,nd 951 logical , intent(in) :: itest 952 integer , intent(in) :: ic,ihyp 953 logical , intent(in) , optional :: est_moy 954 real (kind=dp) , dimension(:,:),allocatable :: triang 955 integer :: stat 956 957 logical :: estmoylocal=.true. 958 959 ! NOTE Olivier : j ai mis est_moy optional et j initlise une variable local 960 ! modlin de jm n a pas besoin d etre modifie, c est donc moin intrusif.... 961 if ( present(est_moy) ) then 962 estmoylocal=est_moy 963 end if 964 965 !cccccccccccccccccccccccc 966 !cccccccccccccccccccccccc 967 ! L'analyse est faite pour le ic �me caract�re 968 ! 969 ! Les premier �l�ment du vecteur des effets � estimer sont la moyenne 970 ! g�n�rale, les effets p�re et les effets m�re 971 ! 972 ! 973 ! initilisation des compteurs 974 ! 975 ! nbfem est le nombre de m�res "estimables" d'apr�s ndmin 976 ! nbef est le nombre d'effets fix�s parasites 977 ! nbco le nombre de covariables 978 ! nbintp le nombre d'effet hierarchis�s dnas les effets qtl p�re 979 ! nbintm idem pour les effets qtl m�re 980 ! nbniv est le nombre total de niveaux d'effets parasites 981 ! nbnivp, celui des efets hierarchis�s dans les effets qtl p�re 982 ! nbnivm, idem pour les qtl m�re 983 ! 984 ! 985 nbfem=0 986 do jm=1,nm 987 if(estime(ic,jm)) then 988 nbfem=nbfem+1 989 indestim(jm)=nbfem 990 end if 991 end do 992 nbef=modele(ic,1) 993 nbtef=nbef 994 nbco=modele(ic,2) 995 nbtco=nbco 996 nbint=modele(ic,3) 997 nbtint=nbint 998 ! 999 ! si itest est � true, il s'agit de retirer un des effet du modele 1000 ! et de recalculer la matrice X'X et les estimabilit� 1001 ! 1002 if(.not.itest) then 1003 meff=0 1004 mcov=0 1005 mint=0 1006 else 1007 if(meff.ne.0) nbef=nbef-1 1008 if(mcov.ne.0) nbco=nbco-1 1009 if(mint.ne.0) nbint=nbint-1 1010 end if 1011 1012 ! 1013 ! construction de la matrice d' incidence xinc 1014 ! de dimension nb de descendants x 1+nb pere + nb mere + 1015 ! nb niveaux effets fix�s + nb c covariables (sous H0) 1016 ! 1017 ! 1018 ! le tableau nivdir donne pour chaque descendant la position 1019 ! dans le vecteur des effets fix�s, des effets exprim�s par ce descendant 1020 ! 1021 ! le tableau covdir donne pour chaque descendant la valeur des covariables 1022 ! retenues dans le mod�le 1023 ! 1024 ! ntniv est le nombre total de niveaux d'effets et covariables 1025 ! 1026 nkd=0 1027 xinc=0.d0 1028 nivdir = 0 1029 do ip=1,np 1030 nm1=nmp(ip)+1 1031 nm2=nmp(ip+1) 1032 do jm=nm1,nm2 1033 nd1=ndm(jm)+1 1034 nd2=ndm(jm+1) 1035 do kd=nd1,nd2 1036 if(presentc(ic,kd)) then 1037 nkd=nkd+1 1038 ! 1039 ! le premier effet fix� est la moyenne g�n�rale 1040 ! 1041 if (estmoylocal) then 1042 nivdir(kd,1)=1 1043 xinc(nkd,1)=1 1044 ntniv=1 1045 nteff=1 1046 else ! pas d esitmation de la moyenne generale 1047 ntniv=0 1048 nteff=0 1049 endif 1050 1051 ! 1052 ! puis, sous H1, les effets QTL 1053 ! 1054 if(ihyp.eq.1) then 1055 ! 1056 ! si des effets fix�s sont en interaction avec l'all�le paternel au qtl 1057 ! il faut estimer les effets qtl pour chacun des niveaux concern�s 1058 ! 1059 ! on commence donc par cr�er un effet composites regroupant l'ensemble de 1060 ! ces effets 1061 ! 1062 nivintp=1 1063 ntlevp=1 1064 ! 07/04/2010 1065 ! correction bug ofi : nbt=3+nbef+nbco n est pas valide 1066 nbt=3+nbtef+nbtco 1067 if(nbint.ne.0) then 1068 ief=0 1069 do jef=1,nbtint 1070 if(mint.ne.jef)then 1071 ief=ief+1 1072 nivintp=nivintp+ntlevp*(nivx(kd,modele(ic,nbt+jef))-1) 1073 ntlevp=ntlevp*nlev(modele(ic,nbt+jef)) 1074 end if 1075 end do 1076 end if 1077 1078 ! 1079 ! puis on cr�e le tableau des niveaux cumul�s et la matrice d'incidence 1080 ! 1081 nivdir(kd,nteff+1)=ntniv+nivintp+ntlevp*(ip-1) 1082 xinc(nkd,nivdir(kd,nteff+1))=prbp(kd) 1083 ntniv=ntniv+ntlevp*np 1084 nteff=nteff+1 1085 ! 1086 ! m�me op�ration pour le qtl transmis par la m�re 1087 ! 1088 if(nbfem.ne.0)then 1089 nivintm=1 1090 ntlevm=1 1091 ! 07/04/2010 1092 ! correction bug ofi : nbt=nbt+nbint n est pas valide 1093 !nbt=nbt+nbint 1094 if(nbint.ne.0) then 1095 ief=0 1096 do jef=1,nbtint 1097 if(nbint.ne.0) then 1098 ief=ief+1 1099 nivintm=nivintm+ntlevm*(nivx(kd,modele(ic,nbt+jef))-1) 1100 ntlevm=ntlevm*nlev(modele(ic,nbt+jef)) 1101 end if 1102 end do 1103 end if 1104 1105 if(estime(ic,jm)) then 1106 km=iam(ic,repfem(jm)) 1107 nivdir(kd,nteff+1)=ntniv+nivintm+ntlevm*(km-1) 1108 xinc(nkd,nivdir(kd,nteff+1))=prbm(kd) 1109 end if 1110 ntniv=ntniv+ntlevm*namest(ic) 1111 nteff=nteff+1 1112 end if 1113 ! 1114 ! fin de H1 1115 ! 1116 end if 1117 ! 1118 ! puis les effets polyg�niques parentaux 1119 ! 1120 if(np.GE.2) then 1121 nivdir(kd,nteff+1)= ntniv+ip 1122 xinc(nkd,nivdir(kd,nteff+1))=1 1123 ntniv=ntniv+np 1124 nteff=nteff+1 1125 endif 1126 if(nbfem.GE.2)then 1127 if(estime(ic,jm)) then 1128 nivdir(kd,nteff+1)=ntniv+indestim(jm) 1129 xinc(nkd,nivdir(kd,nteff+1))=1 1130 end if 1131 ntniv=ntniv+nbfem 1132 nteff=nteff+1 1133 end if 1134 1135 ! 1136 ! autres effets fix�s 1137 ! 1138 nbniv=0 1139 ief=0 1140 do jef=1,nbtef 1141 if(meff.ne.jef) then 1142 ief=ief+1 1143 nivdir(kd,nteff+ief)=ntniv+nbniv+nivx(kd,modele(ic,3+jef)) 1144 xinc(nkd,nivdir(kd,nteff+ief))=1 1145 nbniv=nbniv+nlev(modele(ic,3+jef)) 1146 end if 1147 end do 1148 1149 ntnifix=ntniv+nbniv 1150 ! 07/04/2010 1151 ! correction bug ofi : nbt=3+nbef n est pas valide 1152 nbt=3+nbtef 1153 1154 ! 1155 ! covariables 1156 ! 1157 icov=0 1158 do jcov=1,nbtco 1159 if(mcov.ne.jcov) then 1160 icov=icov+1 1161 covdir(kd,icov)=covar(kd,modele(ic,nbt+jcov)) 1162 xinc(nkd,ntnifix+icov)=covar(kd,modele(ic,nbt+jcov)) 1163 end if 1164 end do 1165 ntniv=ntnifix+nbco 1166 nbt=nbt+nbco 1167 1168 end if 1169 1170 1171 end do 1172 end do 1173 end do 1174 1175 allocate (triang(ntniv,ntniv),STAT=stat) 1176 call check_allocate(stat,'triang [m_qtlmap_analyse_lin_gen]') 1177 1178 ! 1179 ! construction de X'X 1180 ! 1181 xx=0.d0 1182 1183 do kd=1,nkd 1184 do iniv =1,ntniv 1185 do jniv=1,ntniv 1186 xx(iniv,jniv)=xx(iniv,jniv)+xinc(kd,iniv)*xinc(kd,jniv) 1187 end do 1188 end do 1189 end do 1190 1191 ! application de la d�composition de choleski pour d�terminer les contraintes 1192 ! 1193 ! Dans le vecteur vecsol les effets � estimer sont indiqu�s � TRUE et ceux 1194 ! qui ne sont pas estimables � FALSE 1195 ! 1196 ! La matrice triang est la matrice triangulaire sup�rieure M telle que 1197 ! M'M = X'X 1198 ! 1199 do i=1,ntniv 1200 do j=1,ntniv 1201 triang(i,j)=0.d0 1202 end do 1203 end do 1204 1205 do j=1,ntniv 1206 vecsol(j)=.true. 1207 triang(j,j)=xx(j,j) 1208 do k=1,j-1 1209 triang(j,j)=triang(j,j)- triang(j,k)*triang(j,k) 1210 end do 1211 1212 if(triang(j,j).gt.SEUIL_CHO) then 1213 triang(j,j)=dsqrt(triang(j,j)) 1214 do i=j+1,ntniv 1215 triang(i,j)=xx(i,j) 1216 do k=1,j-1 1217 triang(i,j)=triang(i,j)-triang(i,k)*triang(j,k) 1218 end do 1219 triang(i,j)=triang(i,j)/triang(j,j) 1220 end do 1221 1222 else 1223 vecsol(j)=.false. 1224 end if 1225 end do 1226 1227 ! 1228 ! Table de correspondance entre niveaux de d�part et niveaux estimables 1229 ! Au lieu de ntniv effets � estimer, on en a plus que nbnivest 1230 ! corniv(i) est la position, dans la liste des effets estimables 1231 ! du i �me effet initial 1232 ! 1233 1234 nbnivest=0 1235 do ief=1,ntniv 1236 if(vecsol(ief))then 1237 nbnivest=nbnivest+1 1238 corniv(ief)=nbnivest 1239 end if 1240 end do 1241 1242 return 1243 end subroutine contingence_cox
contingence_ldla
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Subroutines ]
NAME
contingence_ldla
DESCRIPTION
INPUTS
ic : index of the trait
ihyp : under this hypothesis
itest : testing nuisances effects
est_moy : estimation of the general mean
chr : chromosome index
n : position tested
details :
est_moy :
option_anal : type analysis * 'LA ', 'LD ', 'LDLA', 'LDJH' *
hsire :
hdam :
OUTPUTS
var_yQy :
NOTES
construction de la matrice de contingence et recherche des contraintes
SOURCE
1271 subroutine contingence_ldla (ic,ihyp,itest,chr,n,var_yQy,details,est_moy,option_anal,hsire,hdam) 1272 ! 1273 ! Tableaux dimensionn�s au nombre d'effets � estimer 1274 integer :: indestim(nm) 1275 ! 1276 ! divers 1277 integer i,j,k,iniv,jniv,ki,kj,kniv,kkd,kd,ief,jef,icov,jcov,chr 1278 integer km,ip,nm1,nm2,nd1,nd2,jm,nbtef,nbtco,nbint,nbtint 1279 integer nkd,nteff,nivintp,nbt,nivintm 1280 integer irank 1281 integer i_gam 1282 logical , intent(in) :: itest 1283 integer , intent(in) :: ic,ihyp,n 1284 logical , intent(in) :: details 1285 logical , intent(in) :: est_moy 1286 character(len=4) , intent(in) :: option_anal 1287 real (kind=dp) var_yQy,ups 1288 real (kind=dp) :: xinc_red(nd,ntnivmax),temp_q(nd,ntnivmax),temp_x(nd,ntnivmax) 1289 real (kind=dp) :: xx_red(ntnivmax,ntnivmax),temp_xx(ntnivmax,ntnivmax) 1290 real (kind=dp) :: mat_q(nd,nd),su 1291 real (kind=dp) , dimension(:,:),allocatable :: triang 1292 integer :: stat 1293 integer i_haplo,geno,j_gam,iall,jp 1294 logical :: estmoylocal,is_la,is_ld,is_ldla,is_ldjh 1295 logical :: hsire,hdam 1296 1297 is_la = (option_anal == 'LA ') 1298 is_ld = (option_anal == 'LD ') 1299 is_ldla = (option_anal == 'LDLA') 1300 is_ldjh = (option_anal == 'LDJH') 1301 1302 1303 estmoylocal=est_moy 1304 1305 !cccccccccccccccccccccccc 1306 !cccccccccccccccccccccccc 1307 ! L'analyse est faite pour le ic �me caract�re 1308 ! 1309 ! Les premier �l�ment du vecteur des effets � estimer sont la moyenne 1310 ! g�n�rale, les effets p�re et les effets m�re 1311 ! 1312 ! 1313 ! initilisation des compteurs 1314 ! 1315 ! nbfem est le nombre de m�res "estimables" d'apr�s ndmin 1316 ! nbef est le nombre d'effets fix�s parasites 1317 ! nbco le nombre de covariables 1318 ! nbintp le nombre d'effet hierarchis�s dnas les effets qtl p�re 1319 ! nbintm idem pour les effets qtl m�re 1320 ! nbniv est le nombre total de niveaux d'effets parasites 1321 ! nbnivp, celui des efets hierarchis�s dans les effets qtl p�re 1322 ! nbnivm, idem pour les qtl m�re 1323 ! 1324 ! 1325 ! print*,'entree dans contingence ldla' 1326 nbfem=0 1327 do jm=1,nm 1328 if(estime(ic,jm)) then 1329 nbfem=nbfem+1 ! correspond a namest 1330 indestim(jm)=nbfem ! correspond a iam ..... 1331 end if 1332 end do 1333 nbef=modele(ic,1) 1334 nbtef=nbef 1335 nbco=modele(ic,2) 1336 nbtco=nbco 1337 nbint=modele(ic,3) 1338 nbtint=nbint 1339 ! 1340 ! si itest est � true, il s'agit de retirer un des effet du modele 1341 ! et de recalculer la matrice X'X et les estimabilit� 1342 ! 1343 if(.not.itest) then 1344 meff=0 1345 mcov=0 1346 mint=0 1347 else 1348 if(meff.ne.0) nbef=nbef-1 1349 if(mcov.ne.0) nbco=nbco-1 1350 if(mint.ne.0) nbint=nbint-1 1351 end if 1352 1353 ! allocate(xinc(nd,size(xx,1))) 1354 1355 ! 1356 ! construction de la matrice d' incidence xinc 1357 ! de dimension nb de descendants x 1+nb pere + nb mere + 1358 ! nb niveaux effets fix�s + nb c covariables (sous H0) 1359 ! 1360 ! 1361 ! le tableau nivdir donne pour chaque descendant la position 1362 ! dans le vecteur des effets fix�s, des effets exprim�s par ce descendant 1363 ! 1364 ! le tableau covdir donne pour chaque descendant la valeur des covariables 1365 ! retenues dans le mod�le 1366 ! 1367 ! ntniv est le nombre total de niveaux d'effets et covariables 1368 ! 1369 nkd=0 1370 xinc=0.d0 1371 nivdir = 0 1372 do ip=1,np 1373 nm1=nmp(ip)+1 1374 nm2=nmp(ip+1) 1375 do jm=nm1,nm2 1376 nd1=ndm(jm)+1 1377 nd2=ndm(jm+1) 1378 do kd=nd1,nd2 1379 if(presentc(ic,kd) ) then 1380 nkd=nkd+1 1381 ! 1382 ! le premier effet fix� est la moyenne g�n�rale 1383 ! 1384 if (estmoylocal) then 1385 nivdir(kd,1)=1 1386 xinc(nkd,1)=1 1387 ntniv=1 1388 nteff=1 1389 else ! pas d estimation de la moyenne generale 1390 ntniv=0 1391 nteff=0 1392 endif 1393 1394 ! 1395 ! puis, sous H1, les effets QTL 1396 ! 1397 if(ihyp.eq.1) then 1398 1399 if(is_ld .or. is_ldla) then 1400 1401 do i_haplo=1,nb_haplo_reduit 1402 pb_haplo(kd,i_haplo)=0.d0 1403 end do 1404 ! 1405 ! haplotype du pere 1406 if(hsire)then 1407 do j=1,2 1408 su = dble(sum(pb_haplo_reduit(num_haplo_pere(ip,j,1:nb_gam_pere(ip,j))))) 1409 if ( su == 0 ) cycle 1410 do i_gam=1,nb_gam_pere(ip,j) 1411 pb_haplo(kd,num_haplo_pere(ip,j,i_gam))= pb_haplo(kd,num_haplo_pere(ip,j,i_gam)) & 1412 +pp_ldla(kd,j) * ( pb_haplo_reduit(num_haplo_pere(ip,j,i_gam)) / su ) 1413 ! print *,'igam:',i_gam,' pb_haplo_reduit:',pb_haplo_reduit(num_haplo_pere(ip,j,i_gam))& 1414 ! ,' pp_ldla:',pp_ldla(kd,j) 1415 end do !i_gam 1416 end do !j 1417 ! 07/04/2010 1418 ! correction bug ofi : nbt=3+nbef+nbco n est pas valide 1419 end if ! opt_model 1420 ! haplotype de la mere si phase connue 1421 if(hdam)then 1422 if(estime(ic,jm)) then 1423 do geno=ngenom(chr,jm)+1,ngenom(chr,jm+1) 1424 do j=1,2 1425 su = dble(sum(pb_haplo_reduit(num_haplo_mere(geno,j,1:nb_gam_mere(geno,j))))) 1426 if ( su == 0 ) cycle 1427 do i_gam=1,nb_gam_mere(geno,j) 1428 pb_haplo(kd,num_haplo_mere(geno,j,i_gam)) = pb_haplo(kd,num_haplo_mere(geno,j,i_gam)) & 1429 +pm_ldla(geno,kd,j) * ( pb_haplo_reduit(num_haplo_mere(geno,j,i_gam)) / su ) 1430 ! print *,'hdam igam:',i_gam,' pb_haplo_reduit:',pb_haplo_reduit(num_haplo_mere(geno,j,i_gam))& 1431 ! ,' pp_ldla:',pm_ldla(geno,kd,j) 1432 end do !i_gam 1433 end do !j 1434 end do !geno 1435 ! haplotype de la mere si phase inconnue 1436 else 1437 ! do j_gam=1,nb_gam(kd) 1438 ! su = dble(sum(pb_haplo_reduit(num_haplo_desc(kd,j_gam,1:nb_gam_desc(kd,j_gam))))) 1439 ! if ( su == 0 ) cycle 1440 ! do i_gam=1,nb_gam_desc(kd,j_gam) 1441 ! pb_haplo(kd,num_haplo_desc(kd,j_gam,i_gam)) = pb_haplo(kd,num_haplo_desc(kd,j_gam,i_gam)) & 1442 ! +prob_gam(kd,j_gam) * (pb_haplo_reduit(num_haplo_desc(kd,j_gam,i_gam)) /su) 1443 ! end do !i_gam 1444 ! end do !j_gam 1445 1446 ! do i_gam=1,nb_gam_desc(kd) 1447 ! pb_haplo(kd,num_haplo_desc(kd,i_gam))= pb_haplo(kd,num_haplo_desc(kd,i_gam))+& 1448 ! pb_haplo_desc(kd,i_gam) 1449 ! end do !i_gam 1450 1451 1452 do i_gam=1,nb_haplo_reduit 1453 pb_haplo(kd,i_gam)= pb_haplo(kd,i_gam)+pb_haplo_desc(kd,i_gam) 1454 end do !i_gam 1455 1456 end if !estime 1457 end if !opt_model (hdam) 1458 ! 1459 ! ligne d'incidence 1460 ! 1461 do i_haplo=1,nb_haplo_reduit 1462 nivdir(kd,nteff+i_haplo)=ntniv+i_haplo 1463 xinc(nkd,nivdir(kd,nteff+i_haplo))=pb_haplo(kd,i_haplo) 1464 ! print *,name_haplo_reduit(i_haplo) 1465 end do 1466 ntniv=ntniv+nb_haplo_reduit 1467 nteff=nteff+nb_haplo_reduit 1468 1469 end if ! option_anal == LD , LDLA 1470 1471 if(option_anal /= 'LD ') then 1472 1473 ! 1474 ! si des effets fix�s sont en interaction avec l'all�le paternel au qtl 1475 ! il faut estimer les effets qtl pour chacun des niveaux concern�s 1476 ! 1477 ! on commence donc par cr�er un effet composites regroupant l'ensemble de 1478 ! ces effets 1479 ! 1480 nivintp=1 1481 ntlevp=1 1482 ! 07/04/2010 1483 ! correction bug ofi : nbt=3+nbef+nbco n est pas valide 1484 nbt=3+nbtef+nbtco 1485 if(nbint.ne.0) then 1486 ief=0 1487 do jef=1,nbtint 1488 if(mint.ne.jef)then 1489 ief=ief+1 1490 nivintp=nivintp+ntlevp*(nivx(kd,modele(ic,nbt+jef))-1) 1491 ntlevp=ntlevp*nlev(modele(ic,nbt+jef)) 1492 end if 1493 end do 1494 end if 1495 1496 ! 1497 ! puis on cr�e le tableau des niveaux cumul�s et la matrice d'incidence 1498 ! 1499 nivdir(kd,nteff+1)=ntniv+nivintp+ntlevp*(ip-1) 1500 xinc(nkd,nivdir(kd,nteff+1))=prbp(kd) 1501 ntniv=ntniv+ntlevp*np 1502 nteff=nteff+1 1503 1504 ! 1505 ! m�me op�ration pour le qtl transmis par la m�re 1506 ! 1507 if(option_anal == 'LDJH') then 1508 if(shared_haplo(kd)/=0) then 1509 jp=int((1+shared_haplo(kd))/2) 1510 iall=shared_haplo(kd)-2*(jp-1) 1511 nivdir(kd,nteff+1)=ntniv+jp 1512 xinc(nkd,nivdir(kd,nteff+1))=xinc(nkd,nivdir(kd,nteff+1))+(-1)**iall 1513 end if 1514 ntniv=ntniv+np 1515 nteff=nteff+1 1516 1517 1518 else 1519 if(nbfem.ne.0)then 1520 nivintm=1 1521 ntlevm=1 1522 ! 07/04/2010 1523 ! correction bug ofi : nbt=nbt+nbint n est pas valide 1524 !nbt=nbt+nbint 1525 if(nbint.ne.0) then 1526 ief=0 1527 do jef=1,nbtint 1528 if(nbint.ne.0) then 1529 ief=ief+1 1530 nivintm=nivintm+ntlevm*(nivx(kd,modele(ic,nbt+jef))-1) 1531 ntlevm=ntlevm*nlev(modele(ic,nbt+jef)) 1532 end if 1533 end do 1534 end if ! nbint ne 0 1535 if(estime(ic,jm)) then 1536 km=iam(ic,repfem(jm)) 1537 nivdir(kd,nteff+1)=ntniv+nivintm+ntlevm*(km-1) 1538 xinc(nkd,nivdir(kd,nteff+1))=prbm(kd) 1539 end if 1540 ntniv=ntniv+ntlevm*namest(ic) 1541 nteff=nteff+1 1542 end if ! nbfem ne 0 1543 end if ! option LDJH 1544 ! 1545 ! fin de l optionLA , LDLA ou LDJH 1546 end if 1547 ! 1548 ! fin de H1 1549 ! 1550 end if 1551 ! 1552 ! puis les effets polyg�niques parentaux 1553 ! 1554 nivdir(kd,nteff+1)= ntniv+ip 1555 xinc(nkd,nivdir(kd,nteff+1))=1 1556 ntniv=ntniv+np 1557 nteff=nteff+1 1558 1559 if(nbfem.ne.0)then 1560 if(estime(ic,jm)) then 1561 nivdir(kd,nteff+1)=ntniv+indestim(jm) 1562 xinc(nkd,nivdir(kd,nteff+1))=1 1563 end if 1564 ntniv=ntniv+nbfem 1565 nteff=nteff+1 1566 end if 1567 1568 ! 1569 ! autres effets fix�s 1570 ! 1571 nbniv=0 1572 ief=0 1573 do jef=1,nbtef 1574 if(meff.ne.jef) then 1575 ief=ief+1 1576 nivdir(kd,nteff+ief)=ntniv+nbniv+nivx(kd,modele(ic,3+jef)) 1577 xinc(nkd,nivdir(kd,nteff+ief))=1 1578 nbniv=nbniv+nlev(modele(ic,3+jef)) 1579 end if 1580 end do 1581 1582 ntnifix=ntniv+nbniv 1583 ! 07/04/2010 1584 ! correction bug ofi : nbt=3+nbef n est pas valide 1585 nbt=3+nbtef 1586 1587 ! 1588 ! covariables 1589 ! 1590 icov=0 1591 do jcov=1,nbtco 1592 if(mcov.ne.jcov) then 1593 icov=icov+1 1594 covdir(kd,icov)=covar(kd,modele(ic,nbt+jcov)) 1595 xinc(nkd,ntnifix+icov)=covar(kd,modele(ic,nbt+jcov)) 1596 end if 1597 end do 1598 ntniv=ntnifix+nbco 1599 nbt=nbt+nbco 1600 1601 end if 1602 1603 end do 1604 end do 1605 end do 1606 1607 allocate (triang(ntniv,ntniv),STAT=stat) 1608 call check_allocate(stat,'triangle [m_qtlmap_analyse_lin_gen]') 1609 1610 ! 1611 ! construction de X'X 1612 ! 1613 xx=0.d0 1614 1615 do kd=1,nkd 1616 do iniv =1,ntniv 1617 do jniv=1,ntniv 1618 xx(iniv,jniv)=xx(iniv,jniv)+xinc(kd,iniv)*xinc(kd,jniv) 1619 end do 1620 end do 1621 end do 1622 1623 ! application de la d�composition de choleski pour d�terminer les contraintes 1624 ! 1625 ! Dans le vecteur vecsol les effets � estimer sont indiqu�s � TRUE et ceux 1626 ! qui ne sont pas estimables � FALSE 1627 ! 1628 ! La matrice triang est la matrice triangulaire sup�rieure M telle que 1629 ! M'M = X'X 1630 ! 1631 do i=1,ntniv 1632 do j=1,ntniv 1633 triang(i,j)=0.d0 1634 end do 1635 end do 1636 1637 do j=1,ntniv 1638 vecsol(j)=.true. 1639 triang(j,j)=xx(j,j) 1640 do k=1,j-1 1641 triang(j,j)=triang(j,j)- triang(j,k)*triang(j,k) 1642 end do 1643 1644 if(triang(j,j).gt.SEUIL_CHO) then 1645 triang(j,j)=dsqrt(triang(j,j)) 1646 do i=j+1,ntniv 1647 triang(i,j)=xx(i,j) 1648 do k=1,j-1 1649 triang(i,j)=triang(i,j)-triang(i,k)*triang(j,k) 1650 end do 1651 triang(i,j)=triang(i,j)/triang(j,j) 1652 end do 1653 1654 else 1655 vecsol(j)=.false. 1656 end if 1657 end do 1658 ! 1659 ! Table de correspondance entre niveaux de d�part et niveaux estimables 1660 ! Au lieu de ntniv effets � estimer, on en a plus que nbnivest 1661 ! corniv(i) est la position, dans la liste des effets estimables 1662 ! du i �me effet initial 1663 ! 1664 1665 nbnivest=0 1666 corniv=0 1667 do ief=1,ntniv 1668 if(vecsol(ief))then 1669 nbnivest=nbnivest+1 1670 corniv(ief)=nbnivest 1671 end if 1672 end do 1673 ! 1674 if(details) then 1675 1676 deallocate (triang) 1677 1678 ! 1679 ! compactage de X'X et xinc 1680 ! 1681 ki=0 1682 do iniv =1,ntniv 1683 if(vecsol(iniv)) then 1684 ki=ki+1 1685 do kd=1,nkd 1686 xinc_red(kd,ki)=xinc(kd,iniv) 1687 end do ! kd 1688 kj=0 1689 do jniv=1,ntniv 1690 if(vecsol(jniv)) then 1691 kj=kj+1 1692 xx_red(ki,kj)=xx(iniv,jniv) 1693 end if 1694 end do ! jniv 1695 end if 1696 end do ! iniv 1697 ! 1698 ! calcul de la variance de la forme quadratique 1699 ! 1700 1701 temp_x=xinc_red(:nkd,:nbnivest) 1702 temp_xx=xx_red(:nbnivest,:nbnivest) 1703 1704 ! 1705 ! inversion de la matrice d'incidence r�duite 1706 ! 1707 ups=1.do-15 1708 call ginv1(temp_xx,nbnivest,size(temp_xx,1),ups,irank) 1709 1710 temp_q=matmul(temp_x,temp_xx) 1711 mat_q(:nbnivest,:nbnivest)=matmul(temp_q,transpose(temp_x)) 1712 1713 var_yQy=0.d0 1714 do iniv=1,nbnivest 1715 var_yQy=var_yQy+1.d0-mat_q(iniv,iniv) 1716 end do 1717 1718 end if ! details 1719 1720 ! deallocate(xinc) 1721 1722 return 1723 end subroutine contingence_ldla
end_analyse_lin_gen
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Subroutines ]
NAME
end_analyse_lin_gen
DESCRIPTION
Deallocation of solution arrays
NOTES
SOURCE
600 subroutine end_analyse_lin_gen 601 deallocate (par0) 602 deallocate (precis0) 603 deallocate (vecsol0) 604 deallocate (par1) 605 deallocate (precis1) 606 deallocate (vecsol1) 607 call end_contingence 608 end subroutine end_analyse_lin_gen
end_contingence
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Subroutines ]
NAME
end_contingence
DESCRIPTION
Deallocation of contingence arrays
NOTES
SOURCE
570 subroutine end_contingence 571 !print *,"** end_contingence ** " 572 deallocate (corniv) 573 deallocate (nivdir) 574 deallocate(xx) 575 deallocate(xxx) 576 deallocate(vecsol) 577 deallocate (covdir) 578 deallocate (prbp) 579 deallocate (prbm) 580 deallocate (ppt) 581 deallocate (pmt) 582 deallocate (pp_ldla) 583 deallocate (pm_ldla) 584 deallocate (pb_haplo) 585 deallocate (prop_haplo_info) 586 deallocate (precis) 587 deallocate(xinc) 588 end subroutine end_contingence
init_analyse_lin_gen
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Subroutines ]
NAME
init_analyse_lin_gen
DESCRIPTION
Initialisation/allocation of solution array
INPUTS
ic : index of trait
nqtl : the hypothesis tested
NOTES
SOURCE
459 subroutine init_analyse_lin_gen(ic,nqtl) 460 integer ,intent(in) :: ic 461 integer ,intent(in) :: nqtl 462 integer :: ntniv,nteff 463 integer :: stat,ngd,npar,chr 464 465 call set_ntnivmax(ic,nqtl,ntnivmax,nteffmax,ntlevp,nbniv) 466 ntlevm=ntlevp 467 npar = np+ntnivmax 468 469 if ( allocated (nmod) ) then 470 if ( nmod(ic)>1) then 471 npar = npar+nmod(ic)-1 472 end if 473 end if 474 475 allocate (par0(npar),STAT=stat) 476 call check_allocate(stat,'par0 [m_qtlmap_analyse_lin_gen]') 477 par0=0.d0 478 allocate (precis0(ntnivmax),STAT=stat) 479 call check_allocate(stat,'precis0 [m_qtlmap_analyse_lin_gen]') 480 precis0=0.d0 481 allocate (vecsol0(ntnivmax),STAT=stat) 482 call check_allocate(stat,'vecsol0 [m_qtlmap_analyse_lin_gen]') 483 484 allocate (par1(npar),STAT=stat) 485 call check_allocate(stat,'par1 [m_qtlmap_analyse_lin_gen]') 486 par1=0.d0 487 allocate (precis1(ntnivmax),STAT=stat) 488 call check_allocate(stat,'precis1 [m_qtlmap_analyse_lin_gen]') 489 precis1=0.d0 490 allocate (vecsol1(ntnivmax),STAT=stat) 491 call check_allocate(stat,'vecsol1 [m_qtlmap_analyse_lin_gen]') 492 493 !Pour garder la compatibilite, on instancie pour le thread courant les 494 !structure de contingence 495 call init_contingence 496 497 end subroutine init_analyse_lin_gen
init_contingence
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Subroutines ]
NAME
init_contingence
DESCRIPTION
Initialisation/allocation of contingence arrays
NOTES
SOURCE
508 subroutine init_contingence 509 integer :: stat,ngd,npar,chr 510 !print *,"** init_contingence ** " 511 allocate (xx(ntnivmax,ntnivmax),STAT=stat) 512 call check_allocate(stat,'xx [init_contingence]') 513 514 allocate (vecsol(ntnivmax),STAT=stat) 515 call check_allocate(stat,'vecsol [init_contingence]') 516 517 allocate (corniv(ntnivmax),STAT=stat) 518 call check_allocate(stat,'corniv [init_contingence]') 519 allocate (nivdir(nd,nteffmax),STAT=stat) 520 call check_allocate(stat,'nivdir [init_contingence]') 521 522 allocate (xxx(ntnivmax,ntnivmax),STAT=stat) 523 call check_allocate(stat,'xxx [init_contingence]') 524 525 allocate (covdir(nd,ncov),STAT=stat) 526 call check_allocate(stat,'covdir [init_contingence]') 527 528 !ngd = 4*nd 529 ! ngd=0 530 ! do chr=1,nchr 531 ! ngd = max(ngd,ngend(chr,ngenom(chr,nm+1)+1)) 532 ! end do 533 !31/08/2010 ! correction bug OFI: on cherche l indice maximum correspondant a ngend 534 ngd = maxval(ngend) 535 536 allocate (pp_ldla(ngd,2),STAT=stat) 537 call check_allocate(stat,'pp_ldla [init_contingence]') 538 !31/08/2010 ! correction bug OFI: size(ngend) => maxval(ngenom) => donne le nombre de genotype possible toute femelles confondues 539 allocate (pm_ldla(maxval(ngenom),ngd,2),STAT=stat) 540 call check_allocate(stat,'pm_ldla [init_contingence]') 541 allocate (pb_haplo(nd,NB_HAPLO_PRIOR),STAT=stat) 542 call check_allocate(stat,'pb_haplo [init_contingence]') 543 allocate (prop_haplo_info(nd),STAT=stat) 544 call check_allocate(stat,'prop_haplo_info [init_contingence]') 545 546 allocate (prbp(nd),STAT=stat) 547 call check_allocate(stat,'prbp [init_contingence]') 548 allocate (prbm(nd),STAT=stat) 549 call check_allocate(stat,'prbm [init_contingence]') 550 allocate (ppt(ngd),STAT=stat) 551 call check_allocate(stat,'ppt [init_contingence]') 552 allocate (pmt(ngd),STAT=stat) 553 call check_allocate(stat,'pmt [init_contingence]') 554 allocate (precis(ntnivmax),STAT=stat) 555 call check_allocate(stat,'precis [init_contingence]') 556 allocate (xinc(nd,ntnivmax),STAT=stat) 557 call check_allocate(stat,'xinc [init_contingence]') 558 end subroutine init_contingence
precision
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Subroutines ]
NAME
precision
DESCRIPTION
INPUTS
xx : incidence matrix
OUTPUTS
precis : vector of precision values
NOTES
SOURCE
2234 subroutine precision(xx,precis) 2235 implicit none 2236 2237 real (kind=dp) ,dimension(ntnivmax,ntnivmax) , intent(in) :: xx 2238 real (kind=dp) ,dimension(ntnivmax), intent(out) :: precis 2239 2240 integer :: iniv,kniv,jniv,lniv,irank 2241 real (kind=dp) :: ups 2242 ! real (kind=dp) ,dimension(ntniv,ntniv) :: xxx 2243 ! 2244 ! R�duction de la matrice d'incidence 2245 ! 2246 xxx=0.d0 2247 2248 kniv=0 2249 do iniv =1,ntniv 2250 if(vecsol(iniv)) then 2251 kniv=kniv+1 2252 lniv=0 2253 do jniv=1,ntniv 2254 if(vecsol(jniv)) then 2255 lniv=lniv+1 2256 xxx(kniv,lniv)=xx(iniv,jniv) 2257 end if 2258 end do 2259 end if 2260 end do 2261 ! 2262 ! inversion de la matrice d'incidence r�duite pour calcul des precisions 2263 ! 2264 2265 ups=1.do-15 2266 call ginv1(xxx,lniv,size(xxx,1),ups,irank) 2267 2268 jniv=0 2269 do iniv=1,ntniv 2270 precis(iniv)=9999.d0 2271 if(vecsol(iniv)) then 2272 jniv=jniv+1 2273 precis(iniv)=xxx(jniv,jniv) 2274 end if 2275 end do 2276 2277 !stop 2278 return 2279 end subroutine precision
prepinc
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Subroutines ]
NAME
prepinc
DESCRIPTION
Prepare the following arrays for the build of the contingence matrix : prbp,prbm,pp_ldla,pm_ldla,pmt,ppt
INPUTS
chr : index of the chromosome
n : the tested position
bcar_icar : the index of trait
option_anal : type analysis * 'LA ', 'LD ', 'LDLA', 'LDJH' *
NOTES
SOURCE
2152 subroutine prepinc(chr,n,bcar_icar,option_anal) 2153 integer , intent(in) :: chr,n 2154 integer , intent(in) :: bcar_icar 2155 character(len=4) ,intent(in) :: option_anal 2156 2157 integer :: ip,nm1,nm2,jm,nd1,nd2,kd,ngeno1,ngeno2,ig,kkd 2158 logical :: is_diff_la=.false. 2159 2160 is_diff_la = (option_anal /= "LA ") 2161 !****************************************************************************** 2162 ! preparation de la matrice d'incidence 2163 ! 2164 ! ************* deb chgt ************* 2165 prbp=0.d0 2166 prbm=0.d0 2167 pp_ldla=0.d0 2168 pm_ldla=0.d0 2169 2170 ! ************* fin chgt ************* 2171 do ip=1,np 2172 nm1=nmp(ip)+1 2173 nm2=nmp(ip+1) 2174 do jm=nm1,nm2 2175 nd1=ndm(jm)+1 2176 nd2=ndm(jm+1) 2177 do kd=nd1,nd2 2178 prbp(kd)=0.d0 2179 prbm(kd)=0.d0 2180 end do 2181 ngeno1=ngenom(chr,jm)+1 2182 ngeno2=ngenom(chr,jm+1) 2183 do ig=ngeno1,ngeno2 2184 nd1=ngend(chr,ig)+1 2185 nd2=ngend(chr,ig+1) 2186 do kd=nd1,nd2 2187 kkd=ndesc(chr,kd) 2188 if(presentc(bcar_icar,kkd)) then 2189 if( (.not. estime(bcar_icar,jm)) .or. ( opt_sib.eq.OPT_SIB_HS) )then 2190 ppt(kd)=-pdd(chr,kd,1,n)+pdd(chr,kd,3,n) 2191 ! ************* deb chgt ************* 2192 if(is_diff_la) then 2193 pp_ldla(kkd,1)=pdd(chr,kd,1,n) 2194 pp_ldla(kkd,2)=pdd(chr,kd,3,n) 2195 end if 2196 ! ************* fin chgt ************* 2197 prbp(kkd)=ppt(kd) 2198 else 2199 ppt(kd)=-pdd(chr,kd,1,n)-pdd(chr,kd,2,n)+pdd(chr,kd,3,n)+pdd(chr,kd,4,n) 2200 prbp(kkd)=prbp(kkd)+probg(chr,ig)*ppt(kd) 2201 pmt(kd)=-pdd(chr,kd,1,n)+pdd(chr,kd,2,n)-pdd(chr,kd,3,n)+pdd(chr,kd,4,n) 2202 prbm(kkd)=prbm(kkd)+probg(chr,ig)*pmt(kd) 2203 ! ************* deb chgt ************* 2204 if(is_diff_la) then 2205 pp_ldla(kkd,1)=pp_ldla(kkd,1)+probg(chr,ig)*(pdd(chr,kd,1,n)+pdd(chr,kd,2,n)) 2206 pp_ldla(kkd,2)=pp_ldla(kkd,2)+probg(chr,ig)*(pdd(chr,kd,3,n)+pdd(chr,kd,4,n)) 2207 pm_ldla(ig,kd,1)=probg(chr,ig)*(pdd(chr,kd,1,n)+pdd(chr,kd,3,n)) 2208 pm_ldla(ig,kd,2)=probg(chr,ig)*(pdd(chr,kd,2,n)+pdd(chr,kd,4,n)) 2209 end if 2210 ! ************* fin chgt ************* 2211 end if 2212 end if 2213 end do 2214 end do 2215 end do 2216 end do 2217 2218 return 2219 end subroutine prepinc
set_filter_optim
[ Top ] [ m_qtlmap_analyse_lin_gen ] [ Subroutines ]
NAME
set_filter_optim
DESCRIPTION
provide a boolean vector that informed which level are expressed according the half sib family
INPUTS
ic : index of the trait
hetero : true => heteroscedastic otherwise homoscedastic
variance_haplo : true => setting the evaluation of a residual variance with the evaluation all haplotypes effects
ntnivmax : maximum number of level
ntniv : number of level of the contingence matrix
vecsol : the boolean vector of the estimability of each effect
xinc : contingence matrix
OUTPUTS
filter_inc : the boolean result array
NOTES
see m_qtlmap_optimization for further details
SOURCE
1747 subroutine set_filter_optim(ic,hetero,variance_haplo,ntnivmax,ntniv,vecsol,xinc,filter_inc) 1748 integer ,intent(in) :: ic,ntnivmax,ntniv 1749 logical ,intent(in) :: hetero,variance_haplo 1750 logical , dimension(ntnivmax) ,intent(in) :: vecsol 1751 real(kind=dp) , dimension(nd,ntnivmax),intent(in) :: xinc 1752 logical , dimension(:,:,:) ,intent(inout) :: filter_inc 1753 1754 integer :: kd1,kd2,ip,jm,ii,i,nstart 1755 1756 !Construction du filtre d optimisation 1757 filter_inc=.true. 1758 if (hetero) then 1759 filter_inc(:,:,1:np)=.false. 1760 nstart = np 1761 else 1762 filter_inc(:,:,1)=.true. 1763 nstart=1 1764 end if 1765 1766 kd1=0 1767 kd2=0 1768 1769 do ip=1,np 1770 if (hetero) then 1771 filter_inc(ip,nmp(ip)+1:nmp(ip+1),ip)=.true. 1772 end if 1773 1774 do jm=nmp(ip)+1,nmp(ip+1) 1775 kd1=kd2+1 1776 kd2=kd2+count(presentc(ic,ndm(jm)+1:ndm(jm+1))) 1777 if (kd1<=kd2) then 1778 ii=0 1779 do i=1,ntniv 1780 if (vecsol(i)) then 1781 ii=ii+1 1782 filter_inc(ip,jm,nstart+ii)=any(xinc(kd1:kd2,i)/=0.d0,dim=1) 1783 end if 1784 end do 1785 else 1786 filter_inc(ip,jm,:)=.false. 1787 end if 1788 end do 1789 end do 1790 1791 1792 1793 1794 end subroutine set_filter_optim