individual.h Source File

00001 //****************************************************
00002 //  April, 1993, University of Illinois
00003 // Copyright (C) 1993, 1994 Tianlin Wang
00004 /* Copyright (C) 1994-2003 Matvec Development Team. 
00005 
00006   This program is free software; you can redistribute it and/or
00007   modify it under the terms of the GNU Library General Public
00008   License as published by the Free Software Foundation; either
00009   version 2 of the License, or (at your option) any later version.
00010   
00011   This program is distributed in the hope that it will be useful,
00012   but WITHOUT ANY WARRANTY; without even the implied warranty of
00013   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00014   Library General Public License for more details.
00015     
00016   You should have received a copy of the GNU Library General Public
00017   License along with this library; if not, write to the Free
00018   Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
00019   MA 02111-1307, USA 
00020 */
00021 
00022 
00023 #ifndef INDIVIDUAL_H
00024 #define INDIVIDUAL_H
00025 
00026 #include <set>
00027 #include <list>
00028 
00029 #include "doublematrix.h"
00030 #include "datanode.h"
00031 #include "genome.h"
00032 #include "geneticdist.h"
00033 #include "gnodederived.h"
00034 #include "dblock.h"
00035 #include "safe_vectors.h"
00036 #include "bg.h"
00037 
00038 namespace matvec {
00039 class Individual;
00040 typedef double (*penetranceType)(const Individual*, const double **);
00041 extern int compare_ind_pt(const void* x, const void* y);
00042 extern int compare_ind_id(const void* a, const void* b);
00043 extern int compare_ind_gid(const void* a, const void* b);
00044 extern int compare_mother_id(const void* x, const void* y);
00045 extern int compare_father_id(const void* x, const void* y);
00046 
00047 class Population;
00048 class Chromosome;
00049 class NuFamily;
00050 
00051 /*!
00052  * \struct RNormal Individual.h
00053  */
00054 struct RNormal
00055 {
00056   double k, a, b, c;
00057   unsigned done;
00058 };
00059 
00060 /*!
00061  * struct UNormal Individual.h
00062  */
00063 struct UNormal
00064 {
00065   double nu,tsq,k;
00066 };
00067 
00068 /*!
00069  * \class Anterior_c
00070  * \brief something here
00071  */
00072 class Anterior_c{
00073 public:
00074   Vector<UNormal> for_gen;
00075 
00076   Anterior_c(void){;};
00077   Anterior_c(unsigned n){for_gen.resize(n);};
00078 };
00079 
00080 /*!
00081  *  \class Posterior_c
00082  *  \brief something here
00083  */
00084 class Posterior_c{
00085 public:
00086   Vector<UNormal> for_gen;
00087   unsigned done;
00088   Posterior_c(void){;};
00089   Posterior_c(unsigned n){for_gen.resize(n); done=0;};
00090 };
00091 
00092 /*!
00093  * \class Mixed_post_vec Individual.h
00094  */
00095 class Mixed_post_vec{
00096  public:
00097   Vector<Vector<UNormal> > postvec;
00098   int done;
00099 };
00100 
00101 extern double global_mu;           //TEMP
00102 extern double global_vg;           //TEMP
00103 extern double global_ve;           //TEMP
00104 
00105 //RLF
00106 /*!
00107  *  \class Individual Individual.h
00108  *  \brief an individual consists of genomes
00109  *
00110  *   \sa Genome
00111  */
00112 class Individual {
00113    friend class Population;
00114    friend class NuFamily;
00115    protected:
00116       unsigned      myid, numchrom, numoffs, numspouse, numtrait;
00117       unsigned      *numoffs_spouse;
00118       Individual    *myfather, *mymother;
00119       Individual    **myoffspring;
00120       Individual    **spouselist;
00121       penetranceType  penetrance_f;
00122 
00123       char          mysex; // changed from an int
00124       int           connect_keeper;
00125 
00126 //BRS
00127   int index_sw, bet_sw, eps_sw; // individual switch index, beta (maternal gamete)  and epsilon (paternal)
00128     int family;
00129       Vector<double>        gprobs;
00130 //BRS
00131       double        inbc;
00132       Genome        genome0, genome1;
00133       Vector<double>*       genotype_counter;
00134 
00135       std::set<Individual *>  offs_tree;
00136       std::list<NuFamily *>   related_family;
00137 
00138       void copyfrom(const Individual& A);
00139 
00140    public:
00141       int           loop_connector,p_origin;
00142       unsigned      group_id;  // an auxilary group id as working variable
00143       int           genotype_id;
00144       double        RBV,MBV, est_GV, true_GV, xbzu_val;
00145       static Population   *population; 
00146       DataNode     *myrecord;
00147       Vector<double>       anterior, *posterior,posterior_iw;
00148       double             anterior_iw;
00149       doubleMatrix       *residual_var;
00150 
00151       Individual(void);
00152       Individual(Population *P);
00153       Individual(const Individual& A);
00154       ~Individual(void){release();}
00155 
00156       const Individual&  operator=(const Individual& A);
00157 
00158       int     terminal() const {return !numoffs;}
00159       int     base() const {return !(myfather || mymother);}
00160       int     isolated() const {return !(myfather || mymother || numoffs);}
00161 
00162       void    reset_id(const unsigned newid) {myid = newid;}
00163       void    release(void);
00164       void    remodel(Population *P);
00165       void    gamete(Chromosome *C, const int n,const double r) const;
00166       void    display() {save(std::cout);}
00167       void    save(std::ostream& out);
00168       void    xbzu(const double x) {xbzu_val = x;}
00169       void    initial_anterior(const double *gfreq,const unsigned tng,
00170                                const int cond=-1);
00171       void    initial_posterior(const unsigned tng,const int cond=-1);
00172       void    pretend_missing(int on,const double *gfreq,const unsigned tng);
00173       void    setPenetrance(penetranceType p){penetrance_f = p;}
00174 
00175       double      get_penetrance(Vector<double> &pen);
00176       double      xbzu(void) const {return xbzu_val;}
00177       double      genotypic_val(void) const;
00178       double      inbcoef(void) const {return inbc;}
00179       double      father_inbcoef(void) const;
00180       double      mother_inbcoef(void) const;
00181 
00182       unsigned    nchrom(void) const {return numchrom;}
00183       unsigned    family_id(void) {return group_id;}
00184       unsigned    marked(void) const {return group_id;}
00185 
00186       unsigned    id(void) const {return myid;}
00187       unsigned    nspouse(void) const {return numspouse;}
00188       unsigned    noffs(void) const {return numoffs;}
00189       const unsigned* noffs_spouse(void) const {return numoffs_spouse;}
00190       char  sex(void) const {return mysex;}
00191 
00192       unsigned    gid(void) const {return group_id;}
00193       unsigned    paternal_chrom_id(const unsigned c) const
00194                        {return genome0.chromosome[c].id();}
00195       unsigned    maternal_chrom_id(const unsigned c) const
00196                        {return genome1.chromosome[c].id();}
00197 
00198       unsigned    father_id(void) const;
00199       unsigned    mother_id(void) const;
00200       unsigned    father_gid(void) const;
00201       unsigned    mother_gid(void) const;
00202 
00203       DataNode*    record(void) const {return myrecord;}
00204       Chromosome*  paternal_chrom(void) {return genome0.chromosome;}
00205       Chromosome*  maternal_chrom(void) {return genome1.chromosome;}
00206       Individual*  father(void) const {return myfather;}
00207       Individual*  mother(void) const {return mymother;}
00208       Individual** offspring(void) const {return myoffspring;}
00209       Individual** spouses(void) const {return spouselist;}
00210 
00211       void   genotype(const unsigned c,const unsigned l,unsigned gtype[]) const;
00212       void   set_genotype(const unsigned c,const unsigned l,const unsigned a0,
00213                           const unsigned a1);
00214 //BRS
00215       void   set_switch(int Nloci);
00216       void   initial_multi_anterior(doubleMatrix& penetrance);
00217       void   initial_multi_posterior(const int cond);
00218       void   collapse_antpost(void);
00219       int    m_anterior_iw;
00220       Dblock m_anterior;
00221       Vector<Dblock> m_posterior;
00222       double m_anterior_scale;
00223       Vector<double> m_posterior_scale;
00224       unsigned n_switches(void);
00225       double get_m_penetrance(doubleMatrix& pen);
00226       void pretend_multi_missing(int on, doubleMatrix& pen);
00227       void put_family(int i){family = i;}
00228       int  get_family(void ){return family;}
00229       void cal_gprobs(Dblock& post_mat_f);
00230       Vector<double> get_gprobs(void){return gprobs;}
00231       double get_m_posterior(int excJ ,Dblock& post_mat);
00232       Vector<Vector<UNormal> > mix_anterior;
00233       Vector<Mixed_post_vec > mix_posterior;
00234       void initial_multi_m_posterior(const unsigned tn_qtl,const int cond );
00235       void initial_multi_m_anterior(const unsigned tn_qtl);
00236       void pretend_multi_m_missing(int on, int tng);
00237       void cal_m_gprobs(int tng);
00238       void collapse_mix_antpost(void);
00239       double get_mix_penetrance(Vector<double> &pen);
00240       static doubleMatrix g_weight;
00241 //BRS
00242 
00243 //MJS
00244       Vector <unsigned> m_gamete, m_gamete1, p_gamete, p_gamete1;
00245       Vector <unsigned> m_founder, p_founder;
00246       static std::vector <std::vector<double> > vec_prob; //SDK was vector
00247       static std::vector <std::vector<double> > vec_cutsetval; //SDK was vector
00248       unsigned p_counter;
00249       unsigned m_counter;
00250 
00251       void put_gametes(Vector <unsigned> m_g,
00252                        Vector <unsigned> p_g);
00253       void  set_founder_alleles(bool);
00254       Vector <unsigned> &get_founder(std::string, bool);
00255       void update_offsprings_founder_alleles(void);
00256       void t0(unsigned SLlocus, unsigned gender_parent);
00257       void t1(unsigned SLlocus);
00258       void t2a(unsigned SLlocus);
00259       void t2b(unsigned SLlocus);
00260       void apply_SL_transition(const unsigned rule, const unsigned lcs);
00261       void apply_SL_cascade(const unsigned rule, const unsigned SLlocus);
00262       double graph_trans_prob(unsigned, unsigned, unsigned);
00263       int marker_index;
00264       double calc_q(void);
00265       void sample_haplotypes(const unsigned);
00266       double q_prob;
00267       doubleMatrix genotype_freq1;
00268       doubleMatrix genotype_freq2;
00269       int ***mode_matrix;      //this is bad stile
00270       int assigned_founder_allele;
00271       void sample_self(int parent_indicator);
00272 //MJS
00273 
00274 //FVP
00275       
00276       Vector <unsigned> id_pdq, id_pdq1; //it stores the pdq index (1,2,3 or 4) for maternal and paternal allele
00277       Vector<unsigned> get_id_pdq(void);
00278       unsigned ord_heter;
00279       unsigned initial_order1, initial_order2;
00280       void search_heteroz(unsigned qtl);
00281       void get_allele_v1(unsigned lcs, Vector<int> allele_vector1);
00282       // HG
00283       void count_haplotype(void);
00284       void display_freq_haplotype(unsigned);
00285       Matrix<unsigned> m_haplotype, p_haplotype;
00286       void map_pdq(int samples, int interval, double r_aq, double r_bq);
00287       void map_pdq(int samples, int interval, BG r_aq, BG r_bq);
00288       double mm_pdq, mp_pdq, pm_pdq, pp_pdq;
00289       BG bgDMPDQ, bgSMPDQ;
00290       Vector<unsigned> m_counter_map, p_counter_map;
00291       void parents(void);
00292       //HG
00293 
00294       //FVP
00295       
00296       ///////6/16/2003 RLF/LRT /////////
00297         SafeSTLVector<unsigned> matHaplotypeCount;
00298         SafeSTLVector<unsigned> patHaplotypeCount;
00299         Vector <unsigned> m_gameteOld, p_gameteOld;
00300         SafeSTLVector<GenotypeNode> genotNodeVector;
00301         SafeSTLVector<AlleleStateNode>  malleleStateNodeVector;
00302         SafeSTLVector<AlleleStateNode>  palleleStateNodeVector;
00303         SafeSTLVector<AlleleOriginNode> malleleOriginNodeVector;
00304         SafeSTLVector<AlleleOriginNode> palleleOriginNodeVector;        
00305         static SafeSTLVector<unsigned> QTLPosVector;
00306         static unsigned numLoci;
00307         static unsigned currentLocus;
00308         static unsigned nQTL;
00309         static double isqrt2pi;
00310         unsigned allAllelesForLocus(unsigned j);
00311         unsigned allGenosForLocus(unsigned j);
00312         double getAllelePenetrance(void);
00313         double getDisAllelePenetrance(void);
00314         double getAllelePenetrance(unsigned forLocus);
00315         double getDisAllelePenetrance(unsigned forLocus);
00316         double getGenoPenetrance(void);
00317         double getDisGenoPenetrance(void);
00318         double getAlleleMTransmissionProb(void);
00319         double getAllelePTransmissionProb(void);
00320         double getAlleleMRTransmissionProb(unsigned forLocus);
00321         double getAllelePRTransmissionProb(unsigned forLocus);  
00322         double getTransitionProb(void);
00323         double getMTransmissionProb(void);
00324         double getPTransmissionProb(void);
00325         double getTransmissionProb(unsigned iA, unsigned mA, unsigned pA,
00326                                    Vector <unsigned> &gamete);
00327         double getMatthiasTransitionProb(void);
00328         double getMatthiasMTransmissionProb(void);
00329         double getMatthiasPTransmissionProb(void);
00330         double getMatthiasTransmissionProb(unsigned iA, unsigned mA, unsigned pA,
00331                                            Vector <unsigned> &gamete,
00332                                            Vector <unsigned> &gameteOld);
00333         int findLeftLocus(Vector <unsigned> &gamete);
00334         int findRightLocus(Vector <unsigned> &gamete);
00335         void setFounderHaplotype(void);
00336         friend ostream& operator<<(ostream& stream, Individual& i);
00337         void displayGenotypes(unsigned locus);
00338         void displayAlleleVectors(unsigned locus);
00339         void setAlleleStateVectors(unsigned locus);
00340         void setAlleleOriginVectors(unsigned locus);
00341         //HG
00342         Vector<Vector<double> > mHaplotypeFreq , pHaplotypeFreq ;
00343         //HG
00344 }; // end of Individual class
00345 
00346 //RLF
00347 extern RNormal* new_n_posterior_vec(const int m);
00348 //RLF
00349 
00350 } // end of namespace matvec
00351 #endif