    001    package calhoun.analysis.crf.features.supporting.phylogenetic;
    002
    003    import org.apache.commons.logging.Log;
    004    import org.apache.commons.logging.LogFactory;
    005
    006    import calhoun.util.Assert;
    007
    008    public class Kimura80Model extends NucleotideSubstitutionModel {
    009            private static final long serialVersionUID = 4075281397219042249L;
    010
    011            private static final Log log = LogFactory.getLog(Kimura80Model.class);
    012
    013            // Model for nucleotide substitution proposed by Kimura in 1980.
    014            // Models the overall rate of nucleotide substitution and the difference
    015            // between rate of transitions and transversions.
    016
    017            double ts;
    018            double tv;
    019
    020            /* Transitions (Ts) are interchanges between pyrimidines (C  T), or between purines (A  G)
    021             * Transversions (Tv) are interchanges beween purines & pyrimidines.
    022             *
    023             * There are twice as many possible transversions, however transitions occur more
    024             * commonly because of the molecular mechanisms involved.
    025             *
    026             * (from http://www.mun.ca/biology/scarr/Transitions_vs_Transversions.htm)
    027             */
    028
    029            public Kimura80Model(double[] parms) {
    030                    setParameters(parms);
    031            }
    032
    033            public Kimura80Model() {
    034                    // Provide some reasonable seed values
    035                    ts = 0.6;
    036                    tv = 0.2;
    037            }
    038
    039            @Override
    040            public String getEvolutionaryModelName() {
    041                    return "Kimura80_nucleotide_evolution_model";
    042            }
    043
    044            @Override
    045            public double[] getParameters() {
    046                    double[] ret = new double[2];
    047                    ret[0] = ts;
    048                    ret[1] = tv;
    049                    return ret;
    050            }
    051
    052
    053            @Override
    054            public void setParameters(double[] parms) {
    055                    Assert.a(parms.length == 2);
    056                    Assert.a(parms[0] > 0);
    057                    Assert.a(parms[1] > 0);
    058
    059                    ts   = parms[0];
    060                    tv = parms[1];
    061
    062                    // Maybe here I should be using a hashing function, but let's first just write it out;
    063                    double[][] X = new double[4][4];
    064                    int A=0, C=1, G=2, T=3;
    065                    X[A][A]=-2*tv-ts;  X[A][C]=tv;       X[A][G]=ts;       X[A][T]=tv;
    066                    X[C][A]=tv;        X[C][C]=-2*tv-ts; X[C][G]=tv;       X[C][T]=ts;
    067                    X[G][A]=ts;        X[G][C]=tv;       X[G][G]=-2*tv-ts; X[G][T]=tv;
    068                    X[T][A]=tv;        X[T][C]=ts;       X[T][G]=tv;       X[T][T]=-2*tv-ts;
    069
    070                    R.assign(X);
    071            }
    072
    073            @Override
    074            public void summarize() {
    075                    log.debug("Kimura80 Nucleotide substitution rate matrix; ts=" + ts + "  tv=" + tv);
    076            }
    077    }