1 /* 2 * modelliemarkov.h 3 * 4 * Created on: 24/05/2016 5 * Author: Michael Woodhams 6 */ 7 8 #ifndef MODELLIEMARKOV_H_ 9 #define MODELLIEMARKOV_H_ 10 11 #include "modelmarkov.h" 12 13 class ModelLieMarkov: public ModelMarkov { 14 public: 15 ModelLieMarkov(string model_name, PhyloTree *tree, string model_params, StateFreqType freq_type, string freq_params); 16 virtual ~ModelLieMarkov(); 17 18 /** 19 this function is served for model testing 20 @param model_name name of the model 21 @param freq_type state frequency type, can be FREQ_USER_DEFINED, FREQ_EQUAL, FREQ_EMPIRICAL, or FREQ_ESTIMATE 22 */ 23 virtual void init(const char *model_name, string model_params, StateFreqType freq, string freq_params); 24 25 /** 26 start structure for checkpointing 27 */ 28 virtual void startCheckpoint(); 29 30 /** 31 save object into the checkpoint 32 */ 33 virtual void saveCheckpoint(); 34 35 /** 36 restore object from the checkpoint 37 */ 38 virtual void restoreCheckpoint(); 39 40 41 /** 42 write information 43 @param out output stream 44 */ 45 virtual void writeInfo(ostream &out); 46 47 static void getLieMarkovModelInfo(string model_name, string &name, string &full_name, int &model_num, int &symmetry, StateFreqType &def_freq); 48 49 static string getModelInfo(string model_name, string &full_name, StateFreqType &def_freq); 50 51 // DO NOT override this function, because 52 // BQM, 2017-05-02: getNDimFreq should return degree of freedom, which is not included in getNDim() 53 // That's why 0 is returned for FREQ_ESTIMATE, num_states-1 for FREQ_EMPIRICAL 54 // virtual int getNDimFreq(); 55 56 // this is redundant, there is already the same function below 57 // bool isTimeReversible(); 58 59 /** 60 @return TRUE if model is time-reversible, FALSE otherwise 61 */ 62 virtual bool isReversible(); 63 64 65 66 static bool validModelName(string model_name); 67 void setBounds(double *lower_bound, double *upper_bound, bool *bound_check); 68 69 /** 70 decompose the rate matrix into eigenvalues and eigenvectors 71 */ 72 virtual void decomposeRateMatrix(); 73 74 /** decompose rate matrix using closed formula derived by Michael Woodhams */ 75 void decomposeRateMatrixClosedForm(); 76 77 /** decompose rate matrix using Eigen library */ 78 virtual void decomposeRateMatrixEigen3lib(); 79 80 /** 81 compute the transition probability matrix. 82 @param time time between two events 83 @param mixture (optional) class for mixture model 84 @param trans_matrix (OUT) the transition matrix between all pairs of states. 85 Assume trans_matrix has size of num_states * num_states. 86 */ 87 virtual void computeTransMatrix(double time, double *trans_matrix, int mixture = 0); 88 // overrides Optimization::restartParameters 89 bool restartParameters(double guess[], int ndim, double lower[], double upper[], bool bound_check[], int iteration); 90 91 protected: 92 /** 93 Model parameters - cached so we know when they change, and thus when 94 recalculations are needed. 95 96 */ 97 double *model_parameters; 98 99 100 double **basis; 101 int symmetry; // RY->0, WS->1, MK->2 102 int model_num; // 0->1.1, etc to 36->12.12 103 void setBasis(); 104 virtual void setRates(); 105 106 /** 107 this function is served for the multi-dimension optimization. It should pack the model parameters 108 into a vector that is index from 1 (NOTE: not from 0) 109 @param variables (OUT) vector of variables, indexed from 1 110 */ 111 virtual void setVariables(double *variables); 112 113 /** 114 this function is served for the multi-dimension optimization. It should assign the model parameters 115 from a vector of variables that is index from 1 (NOTE: not from 0) 116 @param variables vector of variables, indexed from 1 117 @return TRUE if parameters are changed, FALSE otherwise (2015-10-20) 118 */ 119 virtual bool getVariables(double *variables); 120 121 static void parseModelName(string model_name, int* model_num, int* symmetry); 122 /* 123 * Overrides ModelMarkov::getName(). 124 * Avoids appending +FO to name, as this is implied by how LM models 125 * work. 126 * Minh: you might chose to remove this override, if you like "+FO" 127 * to be on LM model names. 128 */ 129 string getName(); 130 131 /* 132 const static double ***BASES; 133 const static int *MODEL_PARAMS; 134 const static string *SYMMETRY; 135 const static string *MODEL_NAMES; 136 const static int NUM_RATES; 137 const static int NUM_LM_MODELS; 138 */ 139 bool validFreqType(); 140 141 }; 142 #endif /* MODELLIEMARKOV_H_ */ 143