Skip to content

Commit

Permalink
Merge pull request #613 from stevenweaver/2.3.1-alpha
Browse files Browse the repository at this point in the history
Merge 2.3.1 into master
  • Loading branch information
stevenweaver authored Aug 28, 2017
2 parents 974948d + 71a90ca commit 93b540c
Show file tree
Hide file tree
Showing 70 changed files with 8,195 additions and 13,048 deletions.
8 changes: 5 additions & 3 deletions .travis.yml
Original file line number Diff line number Diff line change
Expand Up @@ -20,7 +20,10 @@ cache:
- $HOME/cmake-3.3.2-Linux-x86_64/

env:
- MPI=openmpi
- METHOD_TEST=tests/hbltests/libv3/SLAC.wbf
- METHOD_TEST=tests/hbltests/libv3/SLAC-partitioned.wbf
- METHOD_TEST=tests/hbltests/libv3/FEL.wbf
- METHOD_TEST=tests/hbltests/libv3/MEME.wbf

language: c++
compiler:
Expand Down Expand Up @@ -52,5 +55,4 @@ script:
- ./HYPHYMP LIBPATH=`pwd`/res/ tests/hbltests/libv3/math.bf
- ./HYPHYMP LIBPATH=`pwd`/res/ tests/hbltests/libv3/iofunctions.bf
- ./HYPHYMP LIBPATH=`pwd`/res/ tests/hbltests/libv3/utilityfunctions.bf
- ./HYPHYMP LIBPATH=`pwd`/res/ tests/hbltests/libv3/SLAC.wbf
- ./HYPHYMP LIBPATH=`pwd`/res/ tests/hbltests/libv3/SLAC-partitioned.wbf
- ./HYPHYMP LIBPATH=`pwd`/res/ $METHOD_TEST
4 changes: 4 additions & 0 deletions res/TemplateBatchFiles/BUSTED-SRV.bf
Original file line number Diff line number Diff line change
@@ -1,5 +1,9 @@
RequireVersion ("2.1320141020");

console.log("ERROR: THIS BATCHFILE WILL NOT RUN PROPERLY WITH UPDATED TERMS BASE. QUITTING.");
exit();


_BUSTED_timers = {3,1};
busted.taskTimerStart (2);

Expand Down
497 changes: 0 additions & 497 deletions res/TemplateBatchFiles/BUSTED.bf

This file was deleted.

2 changes: 0 additions & 2 deletions res/TemplateBatchFiles/BranchSiteREL.bf
Original file line number Diff line number Diff line change
Expand Up @@ -85,8 +85,6 @@ codon3x4 = BuildCodonFrequencies (nucCF);

tree_info = trees.LoadAnnotatedTopology (1);

//fprintf (stdout, tree_info, "\n");

Model MGL = (MGMatrixLocal, codon3x4, 0);

tree_info_models = tree_info["model_map"];
Expand Down
26 changes: 14 additions & 12 deletions res/TemplateBatchFiles/LHT.bf
Original file line number Diff line number Diff line change
Expand Up @@ -22,7 +22,8 @@ SelectTemplateModel(filteredData);

/*check if there are any partitions.*/

partitionsFound = partitionsFound;

partitionsFound = Rows (DATA_FILE_PARTITION_MATRIX);

fprintf (stdout,"\n\n",partitionsFound," character sets found\n");
if (partitionsFound==0)
Expand Down Expand Up @@ -56,7 +57,7 @@ bestTree = {partitionsFound,1};
bestLk2 = {partitionsFound,1};
globalLk = {ntree,1};

/*GET LIKELIHOOD FOR TREES FOR ALL PARTITIONS AND CHOOSE THE BEST TREE*/
/*GET LIKELIHOOD FOR TREES FOR ALL PARTITIONS AND CHOOSE THE BEST TREE*/
for (part=0; part<partitionsFound; part=part+1)
{
DataSetFilter dsf=CreateFilter (ds,1,DATA_FILE_PARTITION_MATRIX[1][part]);
Expand All @@ -66,14 +67,14 @@ for (part=0; part<partitionsFound; part=part+1)
MBF = PopulateModelMatrix ("filterModelMatrix", baseFreqs);
Model filterModel = (filterModelMatrix, baseFreqs, MBF);
}

fprintf(stdout,"Analysis for partition ",DATA_FILE_PARTITION_MATRIX[0][part]," given trees:\n");
for (i=0; i<ntree; i=i+1)
{
Tree uniqueTree=NEXUS_FILE_TREE_MATRIX[i][1];
LikelihoodFunction lf = (dsf,uniqueTree);
Optimize (mles,lf);
if (part==0)
if (part==0)
{fprintf(stdout,lf);}
listOfLk[part][i]=mles[1][0];
if (i==0)
Expand All @@ -100,7 +101,7 @@ lk1=0; /*H0: constrained same tree for different partition*/
for (i=0; i<ntree; i=i+1)
{
globalLk[i]=0;
for (part=0; part<partitionsFound; part=part+1)
for (part=0; part<partitionsFound; part=part+1)
{
globalLk[i]=globalLk[i]+listOfLk[part][i];
}
Expand All @@ -120,7 +121,7 @@ for (i=0; i<ntree; i=i+1)
}

fprintf (stdout,"\n| Tree n. |");
for (i=0;i<partitionsFound;i=i+1)
for (i=0;i<partitionsFound;i=i+1)
{
fprintf(stdout," ",DATA_FILE_PARTITION_MATRIX[0][i],"\t|");
}
Expand All @@ -143,7 +144,7 @@ fprintf(stdout,"L1= ",Format(lk2,10,10),"\n");
lnlikDelta= 2*(lk2-lk1);
fprintf(stdout,"LRT= ",Format(lnlikDelta,10,10),"\n\n");

/*PARAMETRIC BOOTSTRAP*/
/*PARAMETRIC BOOTSTRAP*/
count=0;

for (replicate=0;replicate<nreplicates;replicate=replicate+1)
Expand All @@ -165,25 +166,26 @@ for (replicate=0;replicate<nreplicates;replicate=replicate+1)

/*simulate constrained different partition*/
DataSet simData = SimulateDataSet(lf);
DataSetFilter simDataf=CreateFilter (simData,1);
DataSetFilter simDataf=CreateFilter (simData,1);
if (FREQUENCY_SENSITIVE)
{
HarvestFrequencies (baseFreqs, simDataf, 1,1,1);
MBF = PopulateModelMatrix ("filterModelMatrix", baseFreqs);
Model filterModel = (filterModelMatrix, baseFreqs, MBF);
}

/*H0: constrained same tree for different partition*/
/*H0: constrained same tree for different partition*/
Tree globalTree=NEXUS_FILE_TREE_MATRIX[theBestTree][1];
LikelihoodFunction lfsim = (simDataf,globalTree);
Optimize (mles,lfsim);
lk1=lk1+mles[1][0];
lk1=lk1+mles[1][0];

/*H1: unconstrained difffeent tree for different partition*/
for (i=0; i<ntree; i=i+1)
{
Tree uniqueTree=NEXUS_FILE_TREE_MATRIX[i][1];
LikelihoodFunction lfsim2 = (simDataf,uniqueTree);
DataSetFilter simDataf2=CreateFilter (simData,1);
LikelihoodFunction lfsim2 = (simDataf2,uniqueTree);
Optimize (mles,lfsim2);
if (i==0) {bestLk=mles[1][0];} else {bestLk=Max(bestLk,mles[1][0]);}
}
Expand All @@ -198,7 +200,7 @@ for (replicate=0;replicate<nreplicates;replicate=replicate+1)
{
fprintf(stdout,"replicate ",replicate,"- 0.00000 *",simLRT,"\n");
}
if (lnlikDelta>simLRT)
if (lnlikDelta>simLRT)
{
count = count+1;
}
Expand Down
26 changes: 16 additions & 10 deletions res/TemplateBatchFiles/ProteinAnalyses/ProteinGTRFit.bf
Original file line number Diff line number Diff line change
Expand Up @@ -15,15 +15,21 @@ utility.SetEnvVariable ("OPTIMIZATION_PRECISION", 0.1); // for testing purposes.
utility.ToggleEnvVariable ("PRODUCE_OPTIMIZATION_LOG", 1); // for testing purposes

io.DisplayAnalysisBanner({
"info": "Fit a general time reversible model to a collection
terms.io.info: "Fit a general time reversible model to a collection
of training protein sequence alignments. Generate substitution
and scoring matrices following the procedures described in Nickle et al 2007",
"version": "0.01",
"reference": "Nickle DC, Heath L, Jensen MA, Gilbert PB, Mullins JI, Kosakovsky Pond SL (2007) HIV-Specific Probabilistic Models of Protein Evolution. PLoS ONE 2(6): e503. doi:10.1371/journal.pone.0000503",
"authors": "Sergei L Kosakovsky Pond and Stephanie J Spielman",
"contact": "{spond,stephanie.spielman}@temple.edu"
terms.io.version: "0.01",
terms.io.reference: "Nickle DC, Heath L, Jensen MA, Gilbert PB, Mullins JI, Kosakovsky Pond SL (2007) HIV-Specific Probabilistic Models of Protein Evolution. PLoS ONE 2(6): e503. doi:10.1371/journal.pone.0000503",
terms.io.authors: "Sergei L Kosakovsky Pond and Stephanie J Spielman",
terms.io.contact: "{spond,stephanie.spielman}@temple.edu"
});


protein_gtr.filename_to_index = terms.data.filename_to_index;
protein_gtr.logl = terms.fit.log_likelihood;
protein_gtr.phase = terms.model.phase;


/********************************************** MENU PROMPTS ********************************************************/
/********************************************************************************************************************/

Expand Down Expand Up @@ -139,7 +145,7 @@ mpi.QueueComplete (protein_gtr.queue);


// Sum of the logL from fitted baseline model across each data set
protein_gtr.baseline_fit_logL = math.Sum (utility.Map (utility.Filter (protein_gtr.analysis_results, "_value_", "_value_/protein_gtr.phase_key"), "_value_", "(_value_[protein_gtr.phase_key])['LogL']"));
protein_gtr.baseline_fit_logL = math.Sum (utility.Map (utility.Filter (protein_gtr.analysis_results, "_value_", "_value_/protein_gtr.phase_key"), "_value_", "(_value_[protein_gtr.phase_key])[terms.fit.log_likelihood]"));

io.ReportProgressMessageMD ("Protein GTR Fitter", " * Initial branch length fit",
"Overall Log(L) = " + protein_gtr.baseline_fit_logL);
Expand All @@ -155,14 +161,14 @@ result_key = "REV-Phase-" + protein_gtr.fit_phase;

if (utility.Has (protein_gtr.analysis_results, result_key, None)) {
io.ReportProgressMessageMD ("Protein GTR Fitter", result_key,
"Loaded cached results for '" + result_key + "'. Log(L) = " + (protein_gtr.analysis_results[result_key])["LogL"] );
"Loaded cached results for '" + result_key + "'. Log(L) = " + (protein_gtr.analysis_results[result_key])[terms.fit.log_likelihood] );
protein_gtr.current_gtr_fit = protein_gtr.analysis_results [result_key];
} else {
protein_gtr.current_gtr_fit = protein_gtr.fitGTRtoFileList (utility.Map (utility.Filter (protein_gtr.analysis_results, "_value_", "_value_/'" + protein_gtr.phase_key + "'"), "_value_", "_value_['" + protein_gtr.phase_key + "']"), protein_gtr.current_gtr_fit, result_key, FALSE);
}

// Record logL
protein_gtr.scores + (protein_gtr.analysis_results[result_key])["LogL"];
protein_gtr.scores + (protein_gtr.analysis_results[result_key])[terms.fit.log_likelihood];

/* Now that the initial GTR fit has been performed, we toggle between a GTR fit and a branch length fit under the estimated GTR parameters */
protein_gtr.shared_EFV = (utility.Values (protein_gtr.current_gtr_fit [terms.efv_estimate]))[0];
Expand Down Expand Up @@ -215,8 +221,8 @@ for (;;) {
for (l1 = 0; l1 < 20; l1 += 1) {
for (l2 = l1 + 1; l2 < 20; l2 += 1) {

previous = (previous_Q[ terms.aminoacidRate (models.protein.alphabet[l1],models.protein.alphabet[l2]) ])[terms.MLE];
current = (current_Q[ terms.aminoacidRate (models.protein.alphabet[l1],models.protein.alphabet[l2]) ])[terms.MLE];
previous = (previous_Q[ terms.aminoacidRate (models.protein.alphabet[l1],models.protein.alphabet[l2]) ])[terms.fit.MLE];
current = (current_Q[ terms.aminoacidRate (models.protein.alphabet[l1],models.protein.alphabet[l2]) ])[terms.fit.MLE];

rmse += (previous - current)^2;
N += 1;
Expand Down
Loading

0 comments on commit 93b540c

Please sign in to comment.