Data Structures | Typedefs | Functions
pll.h File Reference

Data structures for tree and model. More...

Go to the source code of this file.

Data Structures

struct  pllInstanceAttr
 
struct  recompVectors
 Stores the recomputation-state of likelihood vectors. More...
 
struct  ent
 ??? More...
 
struct  hashtable
 ???Hash tables More...
 
struct  stringEnt
 
struct  stringHashtable
 
struct  ratec
 Per-site Rate category entry: likelihood per-site and CAT rate applied ??? More...
 
struct  traversalInfo
 Traversal descriptor entry. More...
 
struct  traversalData
 Traversal descriptor. More...
 
struct  branchInfo
 Branch length information. More...
 
struct  linkageData
 Linkage of partitions. More...
 
struct  linkageList
 
struct  noderec
 Tree node record. More...
 
struct  info
 A brief line. More...
 
struct  bInf
 
struct  iL
 
struct  pInfo
 Partition information structure. More...
 
struct  partitionList
 
struct  checkPointState
 Checkpointing states. More...
 
struct  traversalCounter
 
struct  pllInstance
 Tree topology. More...
 
struct  nniMove
 Stores data related to a NNI move. More...
 
struct  partitionType
 
struct  connectRELL
 
struct  topolRELL
 
struct  topolRELL_LIST
 
struct  conntyp
 Connection within a topology. More...
 
struct  topol
 Single Topology. More...
 
struct  ancestralState
 small helper data structure for printing out/downstream use of marginal ancestral probability vectors. More...
 
struct  bestlist
 List of topologies. More...
 
struct  partitionLengths
 This is used to look up some hard-coded data for each data type. More...
 
struct  pllRollbackInfo
 
struct  pllRearrangeAttr
 Structure holding attributes for searching possible tree rearrangements. More...
 
struct  pllRearrangeInfo
 
struct  pllRearrangeList
 
struct  pllAlignmentData
 Generic structure for storing a multiple sequence alignment. More...
 
struct  threadData
 

Typedefs

typedef int boolean
 
typedef struct ent entry
 
typedef unsigned int hashNumberType
 
typedef struct stringEnt stringEntry
 
typedef struct ratec rateCategorize
 Per-site Rate category entry: likelihood per-site and CAT rate applied ??? More...
 
typedef struct noderec node
 
typedef struct noderecnodeptr
 
typedef struct bInf bestInfo
 
typedef struct iL infoList
 
typedef unsigned int parsimonyNumber
 
typedef struct connectRELLconnptrRELL
 
typedef struct conntyp connect
 Connection within a topology.
 
typedef struct conntypconnptr
 

Functions

double log_approx (double input)
 
double exp_approx (double x)
 
boolean isThisMyPartition (partitionList *pr, int tid, int model)
 Check if a partition is assign to a thread/process. More...
 
void printParallelTimePerRegion (void)
 
void pllFinalizeMPI (void)
 
void pllStartPthreads (pllInstance *tr, partitionList *pr)
 
void pllStopPthreads (pllInstance *tr)
 
int lookupWord (char *s, stringHashtable *h)
 
void pllLockMPI (pllInstance *tr)
 Lock the MPI slave processes prior allocating partitions. More...
 
void pllInitMPI (int *argc, char **argv[])
 Sets up the MPI environment. More...
 
void getDataTypeString (pllInstance *tr, pInfo *partitionInfo, char typeOfData[1024])
 
int countTips (nodeptr p, int numsp)
 
entryinitEntry (void)
 
unsigned int precomputed16_bitcount (unsigned int n, char *bits_in_16bits)
 
double pllGetBranchLength (pllInstance *tr, nodeptr p, int partition_id)
 Get the length of a specific branch. More...
 
void pllSetBranchLength (pllInstance *tr, nodeptr p, int partition_id, double bl)
 Set the length of a specific branch. More...
 
size_t discreteRateCategories (int rateHetModel)
 
const partitionLengthsgetPartitionLengths (pInfo *p)
 
boolean getSmoothFreqs (int dataType)
 
const unsigned int * getBitVector (int dataType)
 
int getUndetermined (int dataType)
 
int getStates (int dataType)
 
char getInverseMeaning (int dataType, unsigned char state)
 
double gettime (void)
 
int gettimeSrand (void)
 
double randum (long *seed)
 
void getxnode (nodeptr p)
 Set the orientation of a node. More...
 
void hookup (nodeptr p, nodeptr q, double *z, int numBranches)
 Connect two nodes and assign branch lengths. More...
 
void hookupFull (nodeptr p, nodeptr q, double *z)
 
void hookupDefault (nodeptr p, nodeptr q)
 
boolean whitechar (int ch)
 
void printLog (pllInstance *tr)
 
void initReversibleGTR (pllInstance *tr, partitionList *pr, int model)
 Initialize GTR. More...
 
double LnGamma (double alpha)
 
double IncompleteGamma (double x, double alpha, double ln_gamma_alpha)
 
double PointNormal (double prob)
 
double PointChi2 (double prob, double v)
 
void makeGammaCats (double alpha, double *gammaRates, int K, boolean useMedian)
 Compute the gamma rates. More...
 
void initModel (pllInstance *tr, double **empiricalFrequencies, partitionList *partitions)
 Initialize the model parameters. More...
 
void resetBranches (pllInstance *tr)
 Reset all branch lengths to default values. More...
 
void modOpt (pllInstance *tr, partitionList *pr, double likelihoodEpsilon)
 
void initializePartitionData (pllInstance *localTree, partitionList *localPartitions)
 
void initMemorySavingAndRecom (pllInstance *tr, partitionList *pr)
 
void makeRandomTree (pllInstance *tr)
 
void nodeRectifier (pllInstance *tr)
 
void makeParsimonyTreeFast (pllInstance *tr, partitionList *pr)
 
void allocateParsimonyDataStructures (pllInstance *tr, partitionList *pr)
 
void freeParsimonyDataStructures (pllInstance *tr, partitionList *pr)
 
void parsimonySPR (nodeptr p, partitionList *pr, pllInstance *tr)
 
FILE * myfopen (const char *path, const char *mode)
 
boolean initrav (pllInstance *tr, partitionList *pr, nodeptr p)
 
void initravPartition (pllInstance *tr, nodeptr p, int model)
 
void update (pllInstance *tr, partitionList *pr, nodeptr p)
 Optimize the length of a specific branch. More...
 
void smooth (pllInstance *tr, partitionList *pr, nodeptr p)
 Branch length optimization of subtree. More...
 
void smoothTree (pllInstance *tr, partitionList *pr, int maxtimes)
 Optimize all branch lenghts of a tree. More...
 
void localSmooth (pllInstance *tr, partitionList *pr, nodeptr p, int maxtimes)
 Optimize the branch length of edges around a specific node. More...
 
boolean localSmoothMulti (pllInstance *tr, nodeptr p, int maxtimes, int model)
 
int pllNniSearch (pllInstance *tr, partitionList *pr, int estimateModel)
 Perform an NNI search. More...
 
int NNI (pllInstance *tr, nodeptr p, int swap)
 Perform an NNI move. More...
 
void smoothRegion (pllInstance *tr, partitionList *pr, nodeptr p, int region)
 Optimize branch lengths of region. More...
 
void regionalSmooth (pllInstance *tr, partitionList *pr, nodeptr p, int maxtimes, int region)
 Wrapper function for optimizing the branch length of a region maxtimes times. More...
 
nodeptr removeNodeBIG (pllInstance *tr, partitionList *pr, nodeptr p, int numBranches)
 Split the tree into two components and optimize new branch length. More...
 
nodeptr removeNodeRestoreBIG (pllInstance *tr, partitionList *pr, nodeptr p)
 Split the tree into two components and recompute likelihood. More...
 
boolean insertBIG (pllInstance *tr, partitionList *pr, nodeptr p, nodeptr q)
 Connect two disconnected tree components. More...
 
boolean insertRestoreBIG (pllInstance *tr, partitionList *pr, nodeptr p, nodeptr q)
 Connect two disconnected tree components without optimizing branch lengths. More...
 
boolean testInsertBIG (pllInstance *tr, partitionList *pr, nodeptr p, nodeptr q)
 Test the.
 
void addTraverseBIG (pllInstance *tr, partitionList *pr, nodeptr p, nodeptr q, int mintrav, int maxtrav)
 Recursively traverse tree and test insertion. More...
 
int rearrangeBIG (pllInstance *tr, partitionList *pr, nodeptr p, int mintrav, int maxtrav)
 Compute the best SPR movement. More...
 
void traversalOrder (nodeptr p, int *count, nodeptr *nodeArray)
 
boolean testInsertRestoreBIG (pllInstance *tr, partitionList *pr, nodeptr p, nodeptr q)
 
void restoreTreeFast (pllInstance *tr, partitionList *pr)
 
void pllTreeEvaluate (pllInstance *tr, partitionList *pr, int maxSmoothIterations)
 Optimize branch lenghts and evaluate likelihood of topology. More...
 
void initTL (topolRELL_LIST *rl, pllInstance *tr, int n)
 Initializes space as large as the tree. More...
 
void freeTL (topolRELL_LIST *rl)
 Deallocate the space associated with this structure. More...
 
void restoreTL (topolRELL_LIST *rl, pllInstance *tr, int n, int numBranches)
 
void resetTL (topolRELL_LIST *rl)
 Reset this structure. More...
 
void saveTL (topolRELL_LIST *rl, pllInstance *tr, int index)
 Save. More...
 
int saveBestTree (bestlist *bt, pllInstance *tr, int numBranches)
 Save the current tree in the bestlist structure. More...
 
int recallBestTree (bestlist *bt, int rank, pllInstance *tr, partitionList *pr)
 Restore the best tree from bestlist structure. More...
 
int initBestTree (bestlist *bt, int newkeep, int numsp)
 Initialize a list of best trees. More...
 
void resetBestTree (bestlist *bt)
 
boolean freeBestTree (bestlist *bt)
 
char * Tree2String (char *treestr, pllInstance *tr, partitionList *pr, nodeptr p, boolean printBranchLengths, boolean printNames, boolean printLikelihood, boolean rellTree, boolean finalPrint, int perGene, boolean branchLabelSupport, boolean printSHSupport)
 
void printTopology (pllInstance *tr, partitionList *pr, boolean printInner)
 
int treeReadLen (FILE *fp, pllInstance *tr, boolean readBranches, boolean readNodeLabels, boolean topologyOnly)
 
void treeReadTopologyString (char *treeString, pllInstance *tr)
 
void getStartingTree (pllInstance *tr)
 
double treeLength (pllInstance *tr, int model)
 
double evaluatePartialGeneric (pllInstance *, partitionList *pr, int i, double ki, int _model)
 
void pllEvaluateGeneric (pllInstance *tr, partitionList *pr, nodeptr p, boolean fullTraversal, boolean getPerSiteLikelihoods)
 Evaluate the log likelihood of the tree topology. More...
 
void pllNewviewGeneric (pllInstance *tr, partitionList *pr, nodeptr p, boolean masked)
 Computes the conditional likelihood vectors of all nodes in the subtree rooted at p. More...
 
void pllNewviewGenericAncestral (pllInstance *tr, partitionList *pr, nodeptr p)
 Computes the conditional likelihood vectors of all nodes in the subtree rooted at p and the marginal ancestral probabilities at node p. More...
 
void newviewAncestralIterative (pllInstance *tr, partitionList *pr)
 A very simple iterative function, we only access the conditional likelihood vector at node p. More...
 
void printAncestralState (nodeptr p, boolean printStates, boolean printProbs, pllInstance *tr, partitionList *pr)
 Prints the ancestral state information for a node p to the terminal. More...
 
void makenewzGeneric (pllInstance *tr, partitionList *pr, nodeptr p, nodeptr q, double *z0, int maxiter, double *result, boolean mask)
 Optimize branch length value(s) of a given branch with the Newton-Raphtson procedure. More...
 
void makenewzGenericDistance (pllInstance *tr, int maxiter, double *z0, double *result, int taxon1, int taxon2)
 
double evaluatePartitionGeneric (pllInstance *tr, nodeptr p, int model)
 
void newviewPartitionGeneric (pllInstance *tr, nodeptr p, int model)
 
double evaluateGenericVector (pllInstance *tr, nodeptr p)
 
void categorizeGeneric (pllInstance *tr, nodeptr p)
 
double makenewzPartitionGeneric (pllInstance *tr, nodeptr p, nodeptr q, double z0, int maxiter, int model)
 
boolean isTip (int number, int maxTips)
 Check whether a node is a tip. More...
 
void computeTraversal (pllInstance *tr, nodeptr p, boolean partialTraversal, int numBranches)
 Compute the traversal descriptor of the subtree rooted at p. More...
 
void allocRecompVectorsInfo (pllInstance *tr)
 Allocates memory for recomputation structure. More...
 
void allocTraversalCounter (pllInstance *tr)
 
boolean getxVector (recompVectors *rvec, int nodenum, int *slot, int mxtips)
 Get a pinned slot slot that holds the likelihood vector for inner node nodenum. More...
 
boolean needsRecomp (boolean recompute, recompVectors *rvec, nodeptr p, int mxtips)
 Checks if the likelihood entries at node p should be updated. More...
 
void unpinNode (recompVectors *v, int nodenum, int mxtips)
 Marks node nodenum as unpinnable. More...
 
void protectNode (recompVectors *rvec, int nodenum, int mxtips)
 Locks node nodenum to force it remains availably in memory. More...
 
void computeTraversalInfoStlen (nodeptr p, int maxTips, recompVectors *rvec, int *count)
 Annotes unoriented tree nodes tr with their subtree size. More...
 
void computeFullTraversalInfoStlen (nodeptr p, int maxTips, recompVectors *rvec)
 Annotes all tree nodes tr with their subtree size. More...
 
void printTraversalInfo (pllInstance *tr)
 
void countTraversal (pllInstance *tr)
 
void pllNewviewIterative (pllInstance *tr, partitionList *pr, int startIndex)
 Compute the conditional likelihood for each entry (node) of the traversal descriptor. More...
 
void pllEvaluateIterative (pllInstance *tr, partitionList *pr, boolean getPerSiteLikelihoods)
 Evaluate the log likelihood of a specific branch of the topology. More...
 
void storeExecuteMaskInTraversalDescriptor (pllInstance *tr, partitionList *pr)
 
void storeValuesInTraversalDescriptor (pllInstance *tr, partitionList *pr, double *value)
 
void makenewzIterative (pllInstance *, partitionList *pr)
 Precompute values (sumtable) from the 2 likelihood vectors of a given branch. More...
 
void execCore (pllInstance *, partitionList *pr, volatile double *dlnLdlz, volatile double *d2lnLdlz2)
 Compute first and second derivatives of the likelihood with respect to a given branch length. More...
 
nodeptr findAnyTip (nodeptr p, int numsp)
 
void putWAG (double *ext_initialRates)
 Hardcoded values for the WAG model. More...
 
unsigned int ** initBitVector (int mxtips, unsigned int *vectorLength)
 
hashtableinitHashTable (unsigned int n)
 
void cleanupHashTable (hashtable *h, int state)
 
double convergenceCriterion (hashtable *h, int mxtips)
 
void freeBitVectors (unsigned int **v, int n)
 
void freeHashTable (hashtable *h)
 
stringHashtableinitStringHashTable (hashNumberType n)
 
void addword (char *s, stringHashtable *h, int nodeNumber)
 
void printBothOpen (const char *format,...)
 
void initRateMatrix (pllInstance *tr, partitionList *pr)
 Initialize the substitution rates matrix. More...
 
void bitVectorInitravSpecial (unsigned int **bitVectors, nodeptr p, int numsp, unsigned int vectorLength, hashtable *h, int treeNumber, int function, branchInfo *bInf, int *countBranches, int treeVectorLength, boolean traverseOnly, boolean computeWRF, int processID)
 
unsigned int bitcount_32_bit (unsigned int i)
 
unsigned int bitcount_64_bit (unsigned long i)
 
void perSiteLogLikelihoods (pllInstance *tr, partitionList *pr, double *logLikelihoods)
 
void updatePerSiteRates (pllInstance *tr, partitionList *pr, boolean scaleRates)
 
void restart (pllInstance *tr, partitionList *pr)
 
double getBranchLength (pllInstance *tr, partitionList *pr, int perGene, nodeptr p)
 
boolean isGap (unsigned int *x, int pos)
 Check whether the position pos in bitvector x is a gap. More...
 
boolean noGap (unsigned int *x, int pos)
 Check whether the position pos in bitvector x is NOT a gap. More...
 
pllNewickTreepllNewickParseString (const char *newick)
 Parse a newick tree string. More...
 
pllNewickTreepllNewickParseFile (const char *filename)
 Parse a newick tree file. More...
 
int pllValidateNewick (pllNewickTree *)
 Validate if a newick tree is a valid phylogenetic tree. More...
 
void pllNewickParseDestroy (pllNewickTree **)
 Deallocate newick parser stack structure. More...
 
int pllNewickUnroot (pllNewickTree *t)
 Convert a binary rooted trree to a binary unrooted tree. More...
 
void pllQueuePartitionsDestroy (pllQueue **partitions)
 Destroy queue structure that contains parsed information from a partition file. More...
 
pllQueuepllPartitionParse (const char *filename)
 Parse a partition (model) file. More...
 
void pllPartitionDump (pllQueue *partitions)
 Dump a parsed partition file in the console. More...
 
void pllAlignmentDataDestroy (pllAlignmentData *)
 Deallocates the memory associated with the alignment data structure. More...
 
int pllAlignmentDataDumpFile (pllAlignmentData *, int, const char *)
 Dump the alignment to a file of format fileFormat. More...
 
void pllAlignmentDataDumpConsole (pllAlignmentData *alignmentData)
 Prints the alignment to the console. More...
 
pllAlignmentDatapllInitAlignmentData (int, int)
 Initialize alignment structure fields. More...
 
pllAlignmentDatapllParseAlignmentFile (int fileType, const char *)
 Parse a file that contains a multiple sequence alignment. More...
 
char * pllReadFile (const char *, long *)
 Read the contents of a file. More...
 
int * pllssort1main (char **x, int n)
 
linkageListinitLinkageList (int *linkList, partitionList *pr)
 
int pllLinkAlphaParameters (char *string, partitionList *pr)
 Link alpha parameters across partitions. More...
 
int pllLinkFrequencies (char *string, partitionList *pr)
 Link base frequency parameters across partitions. More...
 
int pllLinkRates (char *string, partitionList *pr)
 Link Substitution matrices across partitions. More...
 
int pllSetSubstitutionRateMatrixSymmetries (char *string, partitionList *pr, int model)
 Set symmetries among parameters in the Q matrix. More...
 
void pllSetFixedAlpha (double alpha, int model, partitionList *pr, pllInstance *tr)
 Set the alpha parameter of the Gamma model to a fixed value for a partition. More...
 
void pllSetFixedBaseFrequencies (double *f, int length, int model, partitionList *pr, pllInstance *tr)
 Set all base freuqncies to a fixed value for a partition. More...
 
int pllSetOptimizeBaseFrequencies (int model, partitionList *pr, pllInstance *tr)
 Set that the base freuqencies are optimized under ML. More...
 
void pllSetFixedSubstitutionMatrix (double *q, int length, int model, partitionList *pr, pllInstance *tr)
 Set all substitution rates to a fixed value for a specific partition. More...
 
nodeptr pllGetRandomSubtree (pllInstance *)
 Get a random subtree. More...
 
void makeParsimonyTree (pllInstance *tr)
 
void pllPartitionsDestroy (pllInstance *, partitionList **)
 free all data structures associated to a partition More...
 
int pllPartitionsValidate (pllQueue *parts, pllAlignmentData *alignmentData)
 Correspondance check between partitions and alignment. More...
 
partitionListpllPartitionsCommit (pllQueue *parts, pllAlignmentData *alignmentData)
 Constructs the proposed partition scheme. More...
 
void pllAlignmentRemoveDups (pllAlignmentData *alignmentData, partitionList *pl)
 Remove duplicate sites from alignment and update weights vector. More...
 
double ** pllBaseFrequenciesGTR (partitionList *pl, pllAlignmentData *alignmentData)
 Compute the empirical base frequencies for all partitions. More...
 
void pllTreeInitTopologyNewick (pllInstance *, pllNewickTree *, int)
 Initializes the PLL tree topology according to a parsed newick tree. More...
 
int pllLoadAlignment (pllInstance *tr, pllAlignmentData *alignmentData, partitionList *, int)
 Load alignment to the PLL instance. More...
 
void pllEmpiricalFrequenciesDestroy (double ***empiricalFrequencies, int models)
 
void pllTreeInitTopologyRandom (pllInstance *tr, int tips, char **nameList)
 Initialize PLL tree with a random topology. More...
 
void pllTreeInitTopologyForAlignment (pllInstance *tr, pllAlignmentData *alignmentData)
 Initialize a tree that corresponds to a given (already parsed) alignment. More...
 
void pllBaseSubstitute (pllAlignmentData *alignmentData, partitionList *partitions)
 Encode the alignment data to the PLL numerical representation. More...
 
void pllDestroyInstance (pllInstance *)
 Deallocate the PLL instance. More...
 
pllInstancepllCreateInstance (pllInstanceAttr *)
 Create the main instance of PLL. More...
 
int pllInitModel (pllInstance *, partitionList *, pllAlignmentData *)
 Initialize partitions according to model parameters. More...
 
void pllComputeRandomizedStepwiseAdditionParsimonyTree (pllInstance *tr, partitionList *partitions)
 Compute a randomized stepwise addition oder parsimony tree. More...
 
int pllOptimizeModelParameters (pllInstance *tr, partitionList *pr, double likelihoodEpsilon)
 Optimize all free model parameters of the likelihood model. More...
 
pllRearrangeListpllCreateRearrangeList (int max)
 Create a list for storing topology rearrangements. More...
 
void pllDestroyRearrangeList (pllRearrangeList **bestList)
 Deallocator for topology rearrangements list. More...
 
void pllRearrangeSearch (pllInstance *tr, partitionList *pr, int rearrangeType, nodeptr p, int mintrav, int maxtrav, pllRearrangeList *bestList)
 Search for rearrangement topologies. More...
 
void pllRearrangeCommit (pllInstance *tr, partitionList *pr, pllRearrangeInfo *rearr, int saveRollbackInfo)
 Commit a rearrangement move. More...
 
int pllRearrangeRollback (pllInstance *tr, partitionList *pr)
 Rollback the last committed rearrangement move. More...
 
void pllClearRearrangeHistory (pllInstance *tr)
 
void optRateCatPthreads (pllInstance *tr, partitionList *pr, double lower_spacing, double upper_spacing, double *lhs, int n, int tid)
 
void pllMasterBarrier (pllInstance *, partitionList *, int)
 A generic master barrier for executing parallel parts of the code. More...
 
void newviewGTRGAMMAPROT_AVX_LG4 (int tipCase, double *x1, double *x2, double *x3, double *extEV[4], double *tipVector[4], int *ex3, unsigned char *tipX1, unsigned char *tipX2, int n, double *left, double *right, int *wgt, int *scalerIncrement, const boolean useFastScaling)
 
void newviewGTRCAT_AVX_GAPPED_SAVE (int tipCase, double *EV, int *cptr, double *x1_start, double *x2_start, double *x3_start, double *tipVector, int *ex3, unsigned char *tipX1, unsigned char *tipX2, int n, double *left, double *right, int *wgt, int *scalerIncrement, const boolean useFastScaling, unsigned int *x1_gap, unsigned int *x2_gap, unsigned int *x3_gap, double *x1_gapColumn, double *x2_gapColumn, double *x3_gapColumn, const int maxCats)
 
void newviewGTRCATPROT_AVX_GAPPED_SAVE (int tipCase, double *extEV, int *cptr, double *x1, double *x2, double *x3, double *tipVector, int *ex3, unsigned char *tipX1, unsigned char *tipX2, int n, double *left, double *right, int *wgt, int *scalerIncrement, const boolean useFastScaling, unsigned int *x1_gap, unsigned int *x2_gap, unsigned int *x3_gap, double *x1_gapColumn, double *x2_gapColumn, double *x3_gapColumn, const int maxCats)
 
void newviewGTRGAMMA_AVX_GAPPED_SAVE (int tipCase, double *x1_start, double *x2_start, double *x3_start, double *extEV, double *tipVector, int *ex3, unsigned char *tipX1, unsigned char *tipX2, const int n, double *left, double *right, int *wgt, int *scalerIncrement, const boolean useFastScaling, unsigned int *x1_gap, unsigned int *x2_gap, unsigned int *x3_gap, double *x1_gapColumn, double *x2_gapColumn, double *x3_gapColumn)
 
void newviewGTRGAMMAPROT_AVX_GAPPED_SAVE (int tipCase, double *x1_start, double *x2_start, double *x3_start, double *extEV, double *tipVector, int *ex3, unsigned char *tipX1, unsigned char *tipX2, int n, double *left, double *right, int *wgt, int *scalerIncrement, const boolean useFastScaling, unsigned int *x1_gap, unsigned int *x2_gap, unsigned int *x3_gap, double *x1_gapColumn, double *x2_gapColumn, double *x3_gapColumn)
 
void newviewGTRCAT_AVX (int tipCase, double *EV, int *cptr, double *x1_start, double *x2_start, double *x3_start, double *tipVector, int *ex3, unsigned char *tipX1, unsigned char *tipX2, int n, double *left, double *right, int *wgt, int *scalerIncrement, const boolean useFastScaling)
 
void newviewGenericCATPROT_AVX (int tipCase, double *extEV, int *cptr, double *x1, double *x2, double *x3, double *tipVector, int *ex3, unsigned char *tipX1, unsigned char *tipX2, int n, double *left, double *right, int *wgt, int *scalerIncrement, const boolean useFastScaling)
 
void newviewGTRGAMMA_AVX (int tipCase, double *x1_start, double *x2_start, double *x3_start, double *EV, double *tipVector, int *ex3, unsigned char *tipX1, unsigned char *tipX2, const int n, double *left, double *right, int *wgt, int *scalerIncrement, const boolean useFastScaling)
 
void newviewGTRGAMMAPROT_AVX (int tipCase, double *x1, double *x2, double *x3, double *extEV, double *tipVector, int *ex3, unsigned char *tipX1, unsigned char *tipX2, int n, double *left, double *right, int *wgt, int *scalerIncrement, const boolean useFastScaling)
 
void newviewGTRCATPROT_AVX (int tipCase, double *extEV, int *cptr, double *x1, double *x2, double *x3, double *tipVector, int *ex3, unsigned char *tipX1, unsigned char *tipX2, int n, double *left, double *right, int *wgt, int *scalerIncrement, const boolean useFastScaling)
 
int virtual_width (int n)
 
void computeAllAncestralVectors (nodeptr p, pllInstance *tr, partitionList *pr)
 

Detailed Description

Data structures for tree and model.

PLL (version 1.0.0) a software library for phylogenetic inference

Copyright (C) 2013 Tomas Flouri and Alexandros Stamatakis

This program is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version.

This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.

You should have received a copy of the GNU General Public License along with this program. If not, see http://www.gnu.org/licenses/.

For any other enquiries send an Email to Tomas Flouri Tomas.nosp@m..Flo.nosp@m.uri@h.nosp@m.-its.nosp@m..org

When publishing work that uses PLL please cite PLL

ABSTRACT

PLL is a highly optimized, parallelized software library to ease the development of new software tools dealing with phylogenetic inference. Among the functions included in PLL are

DOCUMENTATION

Extensive documentation for using PLL is available online at

            http://www.libpll.org

USAGE

To use PLL,

Author
Tomas Flouri
Fernando Izquierdo-Carrasco
Andre Aberer
Alexandros Stamatakis

Typedef Documentation

typedef struct ratec rateCategorize

Per-site Rate category entry: likelihood per-site and CAT rate applied ???

Function Documentation

void addTraverseBIG ( pllInstance tr,
partitionList pr,
nodeptr  p,
nodeptr  q,
int  mintrav,
int  maxtrav 
)

Recursively traverse tree and test insertion.

Recursively traverses the tree structure starting from node q and tests the insertion of the component specified by leaf p at the edge between q and q->back.

Parameters
trPLL instance
prList of partitions
pLeaf node of one tree component
qEndpoint node of the edge to test the insertion
mintravMinimum radius around q to test the insertion
maxtravMaximum radius around q to test the insertion\

+ Here is the call graph for this function:

+ Here is the caller graph for this function:

void allocRecompVectorsInfo ( pllInstance tr)

Allocates memory for recomputation structure.

Todo:
this should not depend on tr (vectorRecomFraction should be a parameter) PLL_TRUE if recomputation is currently applied
void computeFullTraversalInfoStlen ( nodeptr  p,
int  maxTips,
recompVectors rvec 
)

Annotes all tree nodes tr with their subtree size.

Similar to computeTraversalInfoStlen, but does a full traversal ignoring orientation. The minum cost is defined as the minimum subtree size. In general, the closer a vector is to the tips, the less recomputations are required to re-establish its likelihood entries

Parameters
pPointer to node
maxTipsNumber of tips in the tree
rvecRecomputation info

+ Here is the call graph for this function:

+ Here is the caller graph for this function:

void computeTraversal ( pllInstance tr,
nodeptr  p,
boolean  partialTraversal,
int  numBranches 
)

Compute the traversal descriptor of the subtree rooted at p.

Computes the traversal descriptor of the subtree with root p. By traversal descriptory we essentially mean a preorder traversal of the unrooted topology by rooting it at a node p. If partialTraversal is set to PLL_TRUE then subtrees which are oriented correctly (i.e. if root node r of a subtree has r->x == 1) are not included in the traversal descriptor.

Parameters
trPLL instance
pNode assumed to be the root
partialTraversalIf set to PLL_TRUE, then a partial traversal descriptor is computed.
numBranchesNumber of branches (either per-partition branch or joint branch estimate)

+ Here is the call graph for this function:

+ Here is the caller graph for this function:

void computeTraversalInfoStlen ( nodeptr  p,
int  maxTips,
recompVectors rvec,
int *  count 
)

Annotes unoriented tree nodes tr with their subtree size.

This function recursively updates the subtree size of each inner node.

Note
The subtree size of node p->number is the number of nodes included in the subtree where node record p is the virtual root.
Parameters
pPointer to node
maxTipsNumber of tips in the tree
rvecRecomputation info
countNumber of visited nodes

+ Here is the call graph for this function:

+ Here is the caller graph for this function:

void execCore ( pllInstance tr,
partitionList pr,
volatile double *  _dlnLdlz,
volatile double *  _d2lnLdlz2 
)

Compute first and second derivatives of the likelihood with respect to a given branch length.

Parameters
trlibrary instance
_dlnLdlzFirst derivative dl/dlz
_d2lnLdlz2Second derivative d(dl/dlz)/dlz
Warning
makenewzIterative should have been called to precompute tr->partitionData[model].sumBuffer at the given branch
Note
this function actually computes the first and second derivatives of the likelihood for a given branch stored in tr->coreLZ[model] Note that in the parallel case coreLZ must always be broadcasted together with the traversal descriptor, at least for optimizing branch lengths

+ Here is the caller graph for this function:

void freeTL ( topolRELL_LIST rl)

Deallocate the space associated with this structure.

rl This structure

Todo:
fill the description

+ Here is the caller graph for this function:

void getxnode ( nodeptr  p)

Set the orientation of a node.

Sets the orientation of node p. That means, it will reset the orientation p->next->x and p->next->next->x to 0 and of p->x to 1, meaning that the conditional likelihood vector for that node is oriented on p, i.e. the conditional likelihood vector represents the subtree rooted at p and not any other of the two nodes.

Parameters
pNode which we want to orient

+ Here is the caller graph for this function:

boolean getxVector ( recompVectors rvec,
int  nodenum,
int *  slot,
int  mxtips 
)

Get a pinned slot slot that holds the likelihood vector for inner node nodenum.

If node node nodenum is not pinned to any slot yet, the minimum cost replacement strategy is used.

Parameters
vRecomputation info
nodenumnode id
slotslot id
mxtipsNumber of tips in the tree

+ Here is the caller graph for this function:

void hookup ( nodeptr  p,
nodeptr  q,
double *  z,
int  numBranches 
)

Connect two nodes and assign branch lengths.

Connect the two nodes p and q in each partition i with a branch of length z[i]

Parameters
pNode p
qNode q
numBranchesNumber of partitions

+ Here is the caller graph for this function:

int initBestTree ( bestlist bt,
int  newkeep,
int  numsp 
)

Initialize a list of best trees.

Initialize a list that will contain the best newkeep tree topologies, i.e. the ones that yield the best likelihood. Inside the list initialize space for newkeep + 1 topologies of numsp tips. The additional topology is the starting one

Parameters
btPointer to bestlist to be initialized
newkeepNumber of new topologies to keep
numspNumber of species (tips)
Returns
number of tree topology slots in the list (minus the starting one)
Todo:
Is there a reason that this function is so complex? Many of the checks are unnecessary as the function is called only at two places in the code with newkeep=1 and newkeep=20

+ Here is the call graph for this function:

void initModel ( pllInstance tr,
double **  empiricalFrequencies,
partitionList partitions 
)

Initialize the model parameters.

Initialize the model parameters. Specifically

  • Base frequencies
  • Rate matrix
Parameters
trThe PLL instance
empiricalFrequenciesPointer to the empirical frequencies array
partitionsPointer to the partitions structure
Todo:
What is tr->optimizeRateCategoryInvocations = 1 ?

+ Here is the call graph for this function:

+ Here is the caller graph for this function:

void initRateMatrix ( pllInstance tr,
partitionList pr 
)

Initialize the substitution rates matrix.

Initialize the substitution rates matrices for all partitions

Parameters
trThe PLL instance
prList of partitions
Todo:
Do we need the secondary structure and binary? Will we only use GTR? If yes, we could rename this function to initRateMatrixGTR

+ Here is the call graph for this function:

+ Here is the caller graph for this function:

void initReversibleGTR ( pllInstance tr,
partitionList pr,
int  model 
)

Initialize GTR.

Wrapper function for the decomposition of the substitution rates matrix into eigenvectors and eigenvalues

Parameters
trPLL instance
prList of partitions
modelPartition index

+ Here is the call graph for this function:

+ Here is the caller graph for this function:

void initTL ( topolRELL_LIST rl,
pllInstance tr,
int  n 
)

Initializes space as large as the tree.

Parameters
rlRELL
trPLL instance
nNumber of
Todo:
Don't know what is this used for. Something with RELL?

+ Here is the caller graph for this function:

boolean insertBIG ( pllInstance tr,
partitionList pr,
nodeptr  p,
nodeptr  q 
)

Connect two disconnected tree components.

Connect two disconnected components by specifying an internal edge from one component and a leaf from the other component. The internal edge e is the edge between q and q->back. The leaf is specified by p. Edge e is removed and two new edges are created. The first one is an edge between p->next and q, and the second one is between p->next->next and q->back. The new likelihood vector for node p is computed.

Note
The function makes use of the thoroughInsertion flag
Todo:
What is tr->lzi ? What is thorough insertion? Why do we optimize branch lengths that will be removed? Add explanation
pll.png
The diagram shows in blue colors the new edges that are created and in red the edge that is removed

+ Here is the call graph for this function:

+ Here is the caller graph for this function:

boolean insertRestoreBIG ( pllInstance tr,
partitionList pr,
nodeptr  p,
nodeptr  q 
)

Connect two disconnected tree components without optimizing branch lengths.

Connect two disconnected components by specifying an internal edge from one component and a leaf from the other component. The internal edge e is the edge between q and q->back. The leaf is specified by p. Edge e is removed and two new edges are created. The first one is an edge between p->next and q, and the second one is between p->next->next and q->back. The new likelihood vector for node p is computed.

Note
The function makes use of the thoroughInsertion flag
Todo:
What is the difference between this and insertBIG?

+ Here is the call graph for this function:

boolean isGap ( unsigned int *  x,
int  pos 
)
inline

Check whether the position pos in bitvector x is a gap.

Parameters
xA bitvector represented by unsigned integers
posPosition to check in x if it is set (i.e. it is a gap)
Returns
Returns the value of the bit vector (1 if set, 0 if not)

+ Here is the caller graph for this function:

boolean isThisMyPartition ( partitionList localPr,
int  tid,
int  model 
)

Check if a partition is assign to a thread/process.

Checks whether partition model from partition list localPr is assigned to be processed by process/thread with id tid.

Parameters
localTreeLocal PLL instance
tidThread/Process id
modelPartition number

+ Here is the caller graph for this function:

boolean isTip ( int  number,
int  maxTips 
)

Check whether a node is a tip.

Checks whether the node with number number is a tip.

Parameters
numberNode number to be checked
maxTipsNumber of tips in the tree
Returns
PLL_TRUE if tip, PLL_FALSE otherwise

+ Here is the caller graph for this function:

void localSmooth ( pllInstance tr,
partitionList pr,
nodeptr  p,
int  maxtimes 
)

Optimize the branch length of edges around a specific node.

Optimize maxtimes the branch length of all (3) edges around a given node p of the tree of library instance tr.

Parameters
trThe library instance
pThe node around which to optimize the edges
maxtimesNumber of optimization rounds to perform

+ Here is the call graph for this function:

+ Here is the caller graph for this function:

void makeGammaCats ( double  alpha,
double *  gammaRates,
int  K,
boolean  useMedian 
)

Compute the gamma rates.

Compute the gamma rates

Parameters
alphaAlpha parameter
gammaRatesArray where to store the computed gamma rates
KNumber of categories
useMedianBoolean flag whether to use a median or not
Todo:
Document this more.

+ Here is the caller graph for this function:

void makenewzGeneric ( pllInstance tr,
partitionList pr,
nodeptr  p,
nodeptr  q,
double *  z0,
int  maxiter,
double *  result,
boolean  mask 
)

Optimize branch length value(s) of a given branch with the Newton-Raphtson procedure.

Warning
A given branch may have one or several branch length values (up to PLL_NUM_BRANCHES), usually the later refers to partition-specific branch length values. Thus z0 and result represent collections rather than double values. The number of branch length values is given by tr->numBranches
Parameters
trLibrary instance
pOne node that defines the branch (p->z)
qThe other node side of the branch (usually p->back), but the branch length can be estimated even if p and q are not connected, e.g. before the insertion of a subtree.
z0Initial branch length value(s) for the given branch p->z
maxiterMaximum number of iterations in the Newton-Raphson procedure
resultResulting branch length value(s) for the given branch p->z
maskSpecifies if a mask to track partition convergence (tr->partitionConverged) is being used.
See Also
typical values for maxiter are constants iterations and PLL_NEWZPERCYCLE
Note
Requirement: q->z == p->z

+ Here is the call graph for this function:

+ Here is the caller graph for this function:

void makenewzIterative ( pllInstance tr,
partitionList pr 
)

Precompute values (sumtable) from the 2 likelihood vectors of a given branch.

Warning
These precomputations are stored in tr->partitionData[model].sumBuffer, which is used by function execCore
Parameters
trLibrary instance
Warning
the given branch is implicitly defined in tr by these nodes: pNumber = tr->td[0].ti[0].pNumber; qNumber = tr->td[0].ti[0].qNumber;
Note
This function should be called only once at the very beginning of each Newton-Raphson procedure for optimizing barnch lengths. It initially invokes an iterative newview call to get a consistent pair of vectors at the left and the right end of the branch and thereafter invokes the one-time only precomputation of values (sumtable) that can be re-used in each Newton-Raphson iteration. Once this function has been called we can execute the actual NR procedure

+ Here is the call graph for this function:

+ Here is the caller graph for this function:

boolean needsRecomp ( boolean  recompute,
recompVectors rvec,
nodeptr  p,
int  mxtips 
)

Checks if the likelihood entries at node p should be updated.

A node needs update if one of the following holds:

  1. It is not oriented (p->x == 0)
  2. We are applying recomputations and node p is not currently available in RAM
Parameters
recomputePLL_TRUE if recomputation is currently applied
pNode to check whether it is associated with the likelihood vector
mxtipsNumber of tips in the tree

+ Here is the call graph for this function:

+ Here is the caller graph for this function:

void newviewAncestralIterative ( pllInstance tr,
partitionList pr 
)

A very simple iterative function, we only access the conditional likelihood vector at node p.

+ Here is the call graph for this function:

+ Here is the caller graph for this function:

int NNI ( pllInstance tr,
nodeptr  p,
int  swap 
)

Perform an NNI move.

Modify the tree topology of instance tr by performing an NNI (Neighbour Neighbor Interchange) move at node p. Let q be p->back. If swap is set to PLL_NNI_P_NEXT then the subtrees rooted at p->next->back and q->next->back will be swapped. Otherwise, if swap is set to PLL_NNI_P_NEXTNEXT then the subtrees rooted at p->next->next->back and q->next->back are swapped. For clarity, see the illustration.

Parameters
trPLL instance
pNode to use as origin for performing NNI
swapWhich node to use for the NNI move. PLL_NNI_P_NEXT uses node p->next while PLL_NNI_P_NEXTNEXT uses p->next->next
Returns
In case of success PLL_TRUE, otherwise PLL_FALSE
Todo:
Started error checking here. Instead of checking the errors in the specified way, implement a variadic function where we pass the results of each check and the error code we want to assign if there is at least one negative result
nni.png
In case swap is set to PLL_NNI_P_NEXT then the dashed red edge between p and r is removed and the blue edges are created. If swap is set to PLL_INIT_P_NEXTNEXT then the dashed red edge between p and s is removed and the green edges are created. In both cases the black dashed edge is removed

+ Here is the call graph for this function:

+ Here is the caller graph for this function:

boolean noGap ( unsigned int *  x,
int  pos 
)
inline

Check whether the position pos in bitvector x is NOT a gap.

Parameters
xA bitvector represented by unsigned integers
posPosition to check in x if it is NOT set (i.e. it is NOT a gap)
Returns
Returns the value of the bit vector (1 if set, 0 if not)

+ Here is the caller graph for this function:

void pllAlignmentRemoveDups ( pllAlignmentData alignmentData,
partitionList pl 
)

Remove duplicate sites from alignment and update weights vector.

Removes duplicate sites from the alignment given the partitions list and updates the weight vector of the alignment and the boundaries (upper, lower, width) for each partition.

Parameters
alignmentDataThe multiple sequence alignment
plList of partitions
double** pllBaseFrequenciesGTR ( partitionList pl,
pllAlignmentData alignmentData 
)

Compute the empirical base frequencies for all partitions.

Compute the empirical base frequencies for all partitions in the list pl.

Parameters
plPartition list
alignmentDataMultiple sequence alignment
Returns
A list of pl->numberOfPartitions arrays each of size pl->partitionData[i]->states, where i is the i-th partition

+ Here is the call graph for this function:

+ Here is the caller graph for this function:

void pllBaseSubstitute ( pllAlignmentData alignmentData,
partitionList partitions 
)

Encode the alignment data to the PLL numerical representation.

Transforms the alignment to the PLL internal representation by substituting each base with a specific digit.

Parameters
alignmentDataMultiple sequence alignment
partitionsList of partitions

+ Here is the caller graph for this function:

void pllClearRearrangeHistory ( pllInstance tr)

Clears the rearrangements history from PLL instance

Clears the rearrangements rollback information (history) from the PLL instance tr.

Parameters
trPLL instance

+ Here is the caller graph for this function:

void pllComputeRandomizedStepwiseAdditionParsimonyTree ( pllInstance tr,
partitionList partitions 
)

Compute a randomized stepwise addition oder parsimony tree.

Implements the RAxML randomized stepwise addition order algorithm

Todo:
check functions that are invoked for potential memory leaks!
Parameters
trThe PLL instance
partitionsThe partitions

Implements the RAxML randomized stepwise addition order algorithm

Todo:
check functions that are invoked for potential memory leaks!
Parameters
trThe PLL instance
partitionsThe partitions
pllInstance* pllCreateInstance ( pllInstanceAttr attr)

Create the main instance of PLL.

Create an instance of the phylogenetic likelihood library

Parameters
rateHetModelRate heterogeneity model
fastScalingexplain fastScaling here
saveMemoryexplain saveMemory here
useRecomIf set to PLL_TRUE, enables ancestral state recomputation
Todo:
Document fastScaling, rate heterogeneity and saveMemory and useRecom
Note
Do not set saveMemory to when using useRecom as memory saving techniques are not yet implemented for ancestral state recomputation.
Returns
On success returns an instance to PLL, otherwise NULL

Create an instance of the phylogenetic likelihood library

Parameters
rateHetModelRate heterogeneity model
fastScalingexplain fastScaling here
saveMemoryexplain saveMemory here
useRecomIf set to PLL_TRUE, enables ancestral state recomputation
Todo:
Document fastScaling, rate heterogeneity and saveMemory and useRecom
Note
Do not set saveMemory to when using useRecom as memory saving techniques are not yet implemented for ancestral state recomputation.
Returns
On success returns an instance to PLL, otherwise NULL

+ Here is the call graph for this function:

void pllDestroyInstance ( pllInstance tr)

Deallocate the PLL instance.

Deallocates the library instance and all its elements.

Parameters
trThe PLL instance

+ Here is the call graph for this function:

void pllFinalizeMPI ( void  )

Finalize MPI run

Finalizes MPI run by synchronizing all processes (master + slaves) with a barrier so that all free their allocated resources. Then MPI_Finalize () is called.

Todo:
Similarly as with the pllLockMPI function, this should not be called by the user, but it is called implicitly at the end of pllDestroyInstance. Probably this function should be removed and inline code be placed in pllDestroyInstance.

+ Here is the caller graph for this function:

double pllGetBranchLength ( pllInstance tr,
nodeptr  p,
int  partition_id 
)

Get the length of a specific branch.

Get the length of the branch specified by node p and p->back of partition partition_id. The branch length is decoded from the PLL representation.

Parameters
trPLL instance
pSpecifies one end-point of the branch. The other one is p->back
partition_idSpecifies the partition
Returns
The branch length
nodeptr pllGetRandomSubtree ( pllInstance tr)

Get a random subtree.

Returns the root node of a randomly picked subtree of the tree in PLL instance tr. The picked subtree is guaranteed to have height over 1, that is, the direct descendents of the returned (root) node are not tips.

Parameters
trPLL instance
Returns
The root node of the randomly picked subtree

+ Here is the call graph for this function:

int pllInitModel ( pllInstance tr,
partitionList partitions,
pllAlignmentData alignmentData 
)

Initialize partitions according to model parameters.

Initializes partitions according to model parameters.

Parameters
trThe PLL instance
partitionsList of partitions
alignmentDataThe parsed alignment
Returns
Returns PLL_TRUE in case of success, otherwise PLL_FALSE

+ Here is the call graph for this function:

void pllInitMPI ( int *  argc,
char **  argv[] 
)

Sets up the MPI environment.

Calls the MPI_Init function and makes sure all processes store their process ID and the total number of processes, using a barrier.

Note
this should be the first call that is executed in your main method.
Parameters
argcAddress of argc from main
argvAddress of argv from main
int pllLinkAlphaParameters ( char *  string,
partitionList pr 
)

Link alpha parameters across partitions.

Links alpha paremeters across partitions (GAMMA model of rate heterogeneity)

Parameters
stringstring describing the linkage pattern
prList of partitions
Todo:
test behavior/impact/mem-leaks of this when PSR model is used it shouldn't do any harm, but it would be better to check!

Links alpha paremeters across partitions (GAMMA model of rate heterogeneity)

Parameters
stringstring describing the linkage pattern
prList of partitions
Todo:
test behavior/impact/mem-leaks of this when PSR model is used it shouldn't do any harm, but it would be better to check!
int pllLinkFrequencies ( char *  string,
partitionList pr 
)

Link base frequency parameters across partitions.

Links base frequency paremeters across partitions

Parameters
stringstring describing the linkage pattern
prList of partitions
Todo:
semantics of this function not clear yet: right now this only has an effect when we do a ML estimate of base frequencies when we use empirical or model-defined (protein data) base frequencies, one could maybe average over the per-partition frequencies, but the averages would need to be weighted accodring on the number of patterns per partition

Links base frequency paremeters across partitions

Parameters
stringstring describing the linkage pattern
prList of partitions
Todo:
semantics of this function not clear yet: right now this only has an effect when we do a ML estimate of base frequencies when we use empirical or model-defined (protein data) base frequencies, one could maybe average over the per-partition frequencies, but the averages would need to be weighted accodring on the number of patterns per partition
int pllLinkRates ( char *  string,
partitionList pr 
)

Link Substitution matrices across partitions.

Links substitution matrices (Q matrices) across partitions

Parameters
stringstring describing the linkage pattern
prList of partitions
Todo:
re-think/re-design how this is done for protein models

Links substitution matrices (Q matrices) across partitions

Parameters
stringstring describing the linkage pattern
prList of partitions
Todo:
re-think/re-design how this is done for protein models
void pllLockMPI ( pllInstance tr)

Lock the MPI slave processes prior allocating partitions.

MPI slave processes are locked and wait until the master process has read the number of partitions, which it then broadcasts to slaves, effectively unlocking them. The slave processes will then allocate their own data structures and be locked in the likelihood function.

Parameters
trPLL instance
Todo:
This function should not be called by the user. It is called at pllCreateInstance. Probably this function should be removed and inline code be placed in pllCreateInstance.

+ Here is the call graph for this function:

+ Here is the caller graph for this function:

void pllMasterBarrier ( pllInstance tr,
partitionList pr,
int  jobType 
)

A generic master barrier for executing parallel parts of the code.

A generic master barrier through which the master thread/process controls the work job execution. Through the parameter jobType the master instructs the slaves of what type of work they must conduct.

Parameters
trPLL instance
prList of partitions
jobTypeType of job to be conducted

+ Here is the call graph for this function:

+ Here is the caller graph for this function:

void pllNewviewGeneric ( pllInstance tr,
partitionList pr,
nodeptr  p,
boolean  masked 
)

Computes the conditional likelihood vectors of all nodes in the subtree rooted at p.

Compute the conditional likelihood vectors of all nodes in the subtree rooted at node p. The conditional likelihood vector at node p is recomputed regardless of whether the orientation (i.e. p->x) is correct or not, and, recursuvely, the likelihoods at each node in the subtree as needed and if necessary. In case masked is set to PLL_TRUE, the computation will not take place at partitions for which the conditional likelihood has converged (for example as a reult of previous branch length optimization).

Parameters
trPLL instance
prList of partitions
pRoot of the subtree for which we want to recompute the conditional likelihood vectors
maskedIf set to PLL_TRUE, then likelihood vectors of partitions that are converged are not recomputed.

+ Here is the call graph for this function:

+ Here is the caller graph for this function:

void pllNewviewGenericAncestral ( pllInstance tr,
partitionList pr,
nodeptr  p 
)

Computes the conditional likelihood vectors of all nodes in the subtree rooted at p and the marginal ancestral probabilities at node p.

Compute the conditional likelihood vectors of all nodes in the subtree rooted at node p. The conditional likelihood vector at node p is recomputed regardless of whether the orientation (i.e. p->x) is correct or not, and, recursively, the likelihoods at each node in the subtree as needed and if necessary. In addition, the marginal ancestral probability vector for node p is also computed.

Parameters
trPLL instance
prList of partitions
pNode for which we want to compute the ancestral vector
Note
This function is not implemented with the saveMemory technique.

+ Here is the call graph for this function:

void pllNewviewIterative ( pllInstance tr,
partitionList pr,
int  startIndex 
)

Compute the conditional likelihood for each entry (node) of the traversal descriptor.

Computes the conditional likelihood vectors for each entry (node) in the already computed traversal descriptor, starting from the startIndex entry.

Parameters
trPLL instance
prList of partitions
startIndexFrom which node to start computing the conditional likelihood vectors in the traversal descriptor
Note
This function just iterates over the length of the traversal descriptor and computes the conditional likelihhod arrays in the order given by the descriptor. So in a sense, this function has no clue that there is any tree-like structure in the traversal descriptor, it just operates on an array of structs of given length.

+ Here is the call graph for this function:

+ Here is the caller graph for this function:

int pllNniSearch ( pllInstance tr,
partitionList pr,
int  estimateModel 
)

Perform an NNI search.

Modify the tree topology of instance and model parameters tr by performing a NNI (Neighbour Neighbor Interchange) moves p.

Parameters
trPLL instance
prList of partitions
estimateModelDetermine wheter the model parameters should be optimized
Returns
In case of success PLL_TRUE, otherwise PLL_FALSE

+ Here is the call graph for this function:

int pllOptimizeModelParameters ( pllInstance tr,
partitionList pr,
double  likelihoodEpsilon 
)

Optimize all free model parameters of the likelihood model.

Initializes partitions according to model parameters.

Parameters
trThe PLL instance
prList of partitions
likelihoodEpsilonSpecifies up to which epsilon in likelihood values the iterative routine will be optimizing the parameters

+ Here is the call graph for this function:

void pllPartitionsDestroy ( pllInstance tr,
partitionList **  partitions 
)

free all data structures associated to a partition

frees all data structures allocated for this partition

Parameters
partitionsthe pointer to the partition list
tipsnumber of tips in the tree

+ Here is the call graph for this function:

+ Here is the caller graph for this function:

char* pllReadFile ( const char *  filename,
long *  filesize 
)

Read the contents of a file.

Reads the ile filename and return its content. In addition the size of the file is stored in the input variable filesize. The content of the variable filesize can be anything and will be overwritten.

Parameters
filenameName of the input file
filesizeInput parameter where the size of the file (in bytes) will be stored
Returns
Contents of the file

+ Here is the caller graph for this function:

void pllSetBranchLength ( pllInstance tr,
nodeptr  p,
int  partition_id,
double  bl 
)

Set the length of a specific branch.

Set the length of the branch specified by node p and p->back of partition partition_id. The function encodes the branch length to the PLL representation.

Parameters
trPLL instance
pSpecifies one end-point of the branch. The other one is p->back
partition_idSpecifies the partition
blBranch length
void pllSetFixedAlpha ( double  alpha,
int  model,
partitionList pr,
pllInstance tr 
)

Set the alpha parameter of the Gamma model to a fixed value for a partition.

Sets the alpha parameter of the gamma model of rate heterogeneity to a fixed value and disables the optimization of this parameter

Parameters
alphaalpha value
modelIndex of the partition for which we want to set the alpha value
prList of partitions
trLibrary instance for which we want to fix alpha
Todo:
test if this works with the parallel versions

Sets the alpha parameter of the gamma model of rate heterogeneity to a fixed value and disables the optimization of this parameter

Parameters
alphaalpha value
modelIndex of the partition for which we want to set the alpha value
prList of partitions
trLibrary instance for which we want to fix alpha
Todo:
test if this works with the parallel versions

+ Here is the call graph for this function:

void pllSetFixedBaseFrequencies ( double *  f,
int  length,
int  model,
partitionList pr,
pllInstance tr 
)

Set all base freuqncies to a fixed value for a partition.

Sets all base freuqencies of a partition to fixed values and disables ML optimization of these parameters

Parameters
farray containing the base frequencies
lengthlength of array f, this needs to be as long as the number of states in the model, otherwise an assertion will fail!
modelIndex of the partition for which we want to set the frequencies
prList of partitions
trLibrary instance for which we want to fix the base frequencies
Todo:
test if this works with the parallel versions

Sets all base freuqencies of a partition to fixed values and disables ML optimization of these parameters

Parameters
farray containing the base frequencies
lengthlength of array f, this needs to be as long as the number of states in the model, otherwise an assertion will fail!
modelIndex of the partition for which we want to set the frequencies
prList of partitions
trLibrary instance for which we want to fix the base frequencies
Todo:
test if this works with the parallel versions

+ Here is the call graph for this function:

void pllSetFixedSubstitutionMatrix ( double *  q,
int  length,
int  model,
partitionList pr,
pllInstance tr 
)

Set all substitution rates to a fixed value for a specific partition.

Sets all substitution rates of a partition to fixed values and disables ML optimization of these parameters. It will automatically re-scale the relative rates such that the last rate is 1.0

Parameters
farray containing the substitution rates
lengthlength of array f, this needs to be as long as: (s * s - s) / 2, i.e., the number of upper diagonal entries of the Q matrix
modelIndex of the partition for which we want to set/fix the substitution rates
prList of partitions
trLibrary instance for which we want to fix the substitution rates
Todo:
test if this works with the parallel versions

Sets all substitution rates of a partition to fixed values and disables ML optimization of these parameters. It will automatically re-scale the relative rates such that the last rate is 1.0

Parameters
farray containing the substitution rates
lengthlength of array f, this needs to be as long as: (s * s - s) / 2, i.e., the number of upper diagonal entries of the Q matrix
modelIndex of the partition for which we want to set/fix the substitution rates
prList of partitions
trLibrary instance for which we want to fix the substitution rates
Todo:
test if this works with the parallel versions

+ Here is the call graph for this function:

int pllSetOptimizeBaseFrequencies ( int  model,
partitionList pr,
pllInstance tr 
)

Set that the base freuqencies are optimized under ML.

The base freuqencies for partition model will be optimized under ML

Parameters
modelIndex of the partition for which we want to optimize base frequencies
prList of partitions
trLibrary instance for which we want to fix the base frequencies
Todo:
test if this works with the parallel versions

The base freuqencies for partition model will be optimized under ML

Parameters
modelIndex of the partition for which we want to optimize base frequencies
prList of partitions
trLibrary instance for which we want to fix the base frequencies
Todo:
test if this works with the parallel versions

+ Here is the call graph for this function:

int pllSetSubstitutionRateMatrixSymmetries ( char *  string,
partitionList pr,
int  model 
)

Set symmetries among parameters in the Q matrix.

Allows to link some or all rate parameters in the Q-matrix for obtaining simpler models than GTR

Parameters
stringstring describing the symmetry pattern among the rates in the Q matrix
prList of partitions
modelIndex of the partition for which we want to set the Q matrix symmetries
Todo:
nothing

Allows to link some or all rate parameters in the Q-matrix for obtaining simpler models than GTR

Parameters
stringstring describing the symmetry pattern among the rates in the Q matrix
prList of partitions
modelIndex of the partition for which we want to set the Q matrix symmetries
Todo:
nothing
void pllStartPthreads ( pllInstance tr,
partitionList pr 
)

Start PThreads

Start JOINABLE threads by executing pthread_create. The threads are attached to the pllLikelihoodThread function

Parameters
trPLL instance
prList of partitions
Todo:
This function should never be called by the user. It is called implicitly at pllInitModel. Perhaps we should add a check or inline the code

+ Here is the call graph for this function:

+ Here is the caller graph for this function:

void pllStopPthreads ( pllInstance tr)

Stop PThread

Stop threads by pthread_join

Parameters
trPLL instance
Todo:
This function should never be called by the user. It is implicitly called at pllPartitionsDestroy. We should inline the code

+ Here is the caller graph for this function:

void pllTreeEvaluate ( pllInstance tr,
partitionList pr,
int  maxSmoothIterations 
)

Optimize branch lenghts and evaluate likelihood of topology.

Optimize the branch lengths maxSmoothIterations times and evaluate the likelihood of tree. The resulting likelihood is placed in tr->likelihood

Parameters
trThe PLL instance
prList of partitions
maxSmoothIterationsNumber of times to optimize branch lengths

+ Here is the call graph for this function:

+ Here is the caller graph for this function:

void pllTreeInitTopologyForAlignment ( pllInstance tr,
pllAlignmentData alignmentData 
)

Initialize a tree that corresponds to a given (already parsed) alignment.

Initializes the PLL tree such that it corresponds to the given alignment

Todo:
nothing
Parameters
trThe PLL instance
alignmentDataParsed alignment

Initializes the PLL tree such that it corresponds to the given alignment

Todo:
nothing
Parameters
trThe PLL instance
alignmentDataParsed alignment

+ Here is the call graph for this function:

void pllTreeInitTopologyRandom ( pllInstance tr,
int  tips,
char **  nameList 
)

Initialize PLL tree with a random topology.

Initializes the PLL tree with a randomly created topology

Todo:
Perhaps pass a seed?
Parameters
trThe PLL instance
tipsNumber of tips
nameListA set of tips names representing the taxa labels

*Initializes the PLL tree with a randomly created topology

Todo:
Perhaps pass a seed?
Parameters
trThe PLL instance
tipsNumber of tips
nameListA set of tips names representing the taxa labels

+ Here is the call graph for this function:

void printAncestralState ( nodeptr  p,
boolean  printStates,
boolean  printProbs,
pllInstance tr,
partitionList pr 
)

Prints the ancestral state information for a node p to the terminal.

Prints the ancestral state information for a node p to the terminal. The ancestral state sequence, resp. marginal ancestral state probabilities, is printed depending on whether printStates, resp. printProbs, is set to PLL_TRUE.

Parameters
pThe node for which to print the ancestral state sequence
printStatesIf set to PLL_TRUE then the ancestral state sequence is printed
printProbsIf set to PLL_TRUE then the marginal ancestral state probabilities are printed
trPLL instance
prList of partitions
Note
Here one can see how to store the ancestral probabilities in a dedicated data structure

+ Here is the call graph for this function:

void protectNode ( recompVectors rvec,
int  nodenum,
int  mxtips 
)

Locks node nodenum to force it remains availably in memory.

Warning
If a node is available we dont need to recompute it, but we neet to make sure it is not unpinned while buildding the rest of the traversal descriptor, i.e. unpinnable must be PLL_FALSE at this point, it will automatically be set to PLL_TRUE, after the counter post-order instructions have been executed Omitting this call the traversal will likely still work as long as num_allocated_nodes >> log n, but wrong inner vectors will be used at the wrong moment of pllNewviewIterative, careful!
@param rvec 
  Recomputation info

@param nodenum
  Node id that must remain available in memory 

@param mxtips
  Number of tips in the tree

+ Here is the caller graph for this function:

void putWAG ( double *  ext_initialRates)

Hardcoded values for the WAG model.

Fill the ext_initialRates array with hardcoded substitution rates of the WAG model.

Parameters
ext_initialRatesWhere to place the substitution rates

+ Here is the caller graph for this function:

int rearrangeBIG ( pllInstance tr,
partitionList pr,
nodeptr  p,
int  mintrav,
int  maxtrav 
)

Compute the best SPR movement.

Compute all SPR moves starting from p in the space defined by mintrav and maxtrav and store the best in the tr structure.

Parameters
trPLL instancve
prList of partitions
pNode from which to start the SPR moves testing
mintravMinimum distance from p where to start testing SPRs
maxtravMaximum distance from p where to test SPRs
Returns
0,1 or PLL_BADREAR
Todo:
fix the return value

q1->tip

q2->tip

+ Here is the call graph for this function:

+ Here is the caller graph for this function:

int recallBestTree ( bestlist bt,
int  rank,
pllInstance tr,
partitionList pr 
)

Restore the best tree from bestlist structure.

Restore the rank-th best tree from the bestlist structure bt.

Parameters
btThe bestlist structure containing the stored best trees
rankThe rank (by score) of the tree we want to retrieve
trPLL instance
prList of partitions
Returns
Index (rank) of restored topology in bestlist
void regionalSmooth ( pllInstance tr,
partitionList pr,
nodeptr  p,
int  maxtimes,
int  region 
)

Wrapper function for optimizing the branch length of a region maxtimes times.

Optimize the branch lengths of a specific region maxtimes times. The branch optimization starts at a given node p and is carried out in all nodes with distance upto region from p.

Parameters
trThe library instance.
pNode to start branch optimization from.
maxtimesNumber of times to perform branch optimization.
regionThe allwed node distance from for which to still perform branch optimization.
Todo:
In the previous version (before the model-sep merge) the loops were controlled by tr->numBranches, and now they are controlled by a constant PLL_NUM_BRANCHES. What is right?

+ Here is the call graph for this function:

nodeptr removeNodeBIG ( pllInstance tr,
partitionList pr,
nodeptr  p,
int  numBranches 
)

Split the tree into two components and optimize new branch length.

Split the tree into two components. The disconnection point is node p. First, a branch length is computed for the newly created branch between nodes p->next->back and p->next->next->back and then the two nodes are connected (hookup). Disconnection is done by setting p->next->next->back and p->next->back to NULL.

Parameters
trThe library instance
pThe node at which the tree should be decomposed into two components.
numBranchesNumber of branches per partition
Returns
Node from the disconnected component
Todo:
Why do we return this node?
removeBIG.png
The diagram shows in blue color the new edge that is created and in red the edges that are removed

+ Here is the call graph for this function:

+ Here is the caller graph for this function:

nodeptr removeNodeRestoreBIG ( pllInstance tr,
partitionList pr,
nodeptr  p 
)

Split the tree into two components and recompute likelihood.

Split the tree into two component. The disconnection point is node p. Set the branch length of the new node between p->next->back and p->next->next->back to tr->currentZQR and then decompose the tree into two components by setting p->next->back and p->next->next->back to NULL.

Parameters
trThe library instance
pThe node at which the tree should be decomposed into two components.
Returns
q the node after p
Todo:
Why do we return this node? Why do we set to tr->currentZQR and not compute new optimized length? What is tr->currentZQR?

+ Here is the call graph for this function:

void resetBranches ( pllInstance tr)

Reset all branch lengths to default values.

Reset all branch lengths in the tree instance to default values (PLL_DEFAULTZ)

Parameters
trPLL instance

+ Here is the caller graph for this function:

void resetTL ( topolRELL_LIST rl)

Reset this structure.

Reset the likelihoods in this structure

Parameters
rlThis structure
Todo:
Complete this
int saveBestTree ( bestlist bt,
pllInstance tr,
int  numBranches 
)

Save the current tree in the bestlist structure.

Save the current tree topology in bestlist structure bt.

Parameters
trThe PLL instance
btThe bestlist structure
numBranchesNumber of branches u
Returns
it is never used
Todo:
What to do with the return value? Should we simplify the code?

+ Here is the caller graph for this function:

void saveTL ( topolRELL_LIST rl,
pllInstance tr,
int  index 
)

Save.

Save this topology?

Todo:
Complete this

+ Here is the caller graph for this function:

void smooth ( pllInstance tr,
partitionList pr,
nodeptr  p 
)

Branch length optimization of subtree.

Optimize the length of branch connected by p and p->back, and the lengths of all branches in the subtrees rooted at p->next and p->next->next

Parameters
trThe library instance
prPartition list
pEndpoint of branches to be optimized

+ Here is the call graph for this function:

+ Here is the caller graph for this function:

void smoothRegion ( pllInstance tr,
partitionList pr,
nodeptr  p,
int  region 
)

Optimize branch lengths of region.

Optimize the branch lenghts of only a specific region. The branch optimization starts at a node p and is carried out in all nodes with distance upto region edges from p.

Parameters
trThe library instance.
pNode to start branch optimization from.
regionThe allowed node distance from for which to still perform branch optimization.

+ Here is the call graph for this function:

+ Here is the caller graph for this function:

void smoothTree ( pllInstance tr,
partitionList pr,
int  maxtimes 
)

Optimize all branch lenghts of a tree.

Perform maxtimes rounds of branch length optimization by running smooth() on all neighbour nodes of node tr->start.

Parameters
trThe library instance
maxtimesNumber of optimization rounds to perform

+ Here is the call graph for this function:

+ Here is the caller graph for this function:

void unpinNode ( recompVectors v,
int  nodenum,
int  mxtips 
)

Marks node nodenum as unpinnable.

The slot holding the node nodenum is added to the pool of slot candidates that can be overwritten.

Parameters
vRecomputation info
nodenumnode id
mxtipsNumber of tips in the tree

+ Here is the caller graph for this function:

void update ( pllInstance tr,
partitionList pr,
nodeptr  p 
)

Optimize the length of a specific branch.

Optimize the length of the branch connecting p and p->back for each partition (tr->numBranches) in library instance tr.

Parameters
trThe library instance
prPartition list
pEndpoints of branch to be optimized

+ Here is the call graph for this function:

+ Here is the caller graph for this function: