Generate MrBayes block from IQ-TREE on CIPRES

[eosthusiast]

Is there some way to generate lset and prset commands in a MrBayes data block for each partition that was generated using IQ-TREE's built-in PartionFinder3 on CIPRES? I'm basically looking for something that has the functionality of PhyloSuite or partitNex2MrBayes, both of which are still running PartitionFinder2 and are very slow off the cluster.

Specifically, I used these commands in IQ-TREE:
iqtree2 -nt 6 -bnni -s infile.txt -p partition_file.txt -st DNA -mset mrbayes -merit AICC -m MFP --sprrad 6 -B 1000

This was the output in the best_scheme.nex file:

#nexus
begin sets;
charset ITS = 1-795;
charset LSU = 796-1753;
charset RBP1-1 = 1754-1946;
charset R1Intron1 = 1947-2075;
charset RBP1-2 = 2076-2823;
charset R1Intron2 = 2824-2886;
charset RBP1-3 = 2887-3503;
charset RPB2-1 = 3504-3596;
charset RIntron1 = 3597-3661;
charset RPB2-2 = 3662-4171;
charset RIntron2 = 4172-4248;
charset RPB2-3 = 4249-5247;
charset Tef-1 = 5248-5451;
charset TIntron1 = 5452-5518;
charset Tef-2 = 5519-5655;
charset TIntron2 = 5656-5721;
charset Tef-3 = 5722-5864;
charpartition mymodels =
GTR+F+I+G4: ITS,
GTR+F+I+G4: LSU,
SYM+I+R2: RBP1-1,
GTR+F+G4: R1Intron1,
GTR+F+I+G4: RBP1-2,
K2P+I: R1Intron2,
SYM+I+G4: RBP1-3,
K2P+I+G4: RPB2-1,
GTR+F+R2: RIntron1,
SYM+I+G4: RPB2-2,
GTR+F+R2: RIntron2,
SYM+I+R3: RPB2-3,
GTR+F+I+R3: Tef-1,
HKY+F+I+G4: TIntron1,
GTR+F+I+G4: Tef-2,
HKY+F+I+G4: TIntron2,
GTR+F+I+G4: Tef-3;
end;

Thanks in advance!

[protokute]

This would be extremely helpful. Looking forward to any improvement on this.

[roblanf]

@eosthusiast there are two answers here.

The short answer is 'no', there is no way to output a MrBayes block from IQ-TREE.

The slightly longer answer, and maybe more useful one, is that IQ-TREE does output something in nexus format that you could use as a good starting point: the best_model.nex file, and one could borrow from the PartitionFinder code to implement a conversion script:

Here's an example of the best_model.nex file:

#nexus
begin sets;
  charset dloop = , 1-1114;
  charset a61_a81 = , 1115-1792\3  1793-1951\3;
  charset a62 = , 1116-1792\3;
  charset a63_a83_c13_c23_c33_cytb3_n13_n23_n33_n43_n4l3_n53 = , 1117-1792\3  1795-1951\3  1954-3490\3  3493-4171\3  4174-4951\3  4954-6085\3  6088-7039\3  7042-8080\3  8083-8425\3  8428-9802\3  9805-10090\3  10093-11899\3;
  charset a82_cytb1_n11_n21_n51_12s_loops_16s_loops = , 1794-1951\3  4952-6085\3  6086-7039\3  7040-8080\3  10091-11899\3  12422-12903  12904-13857;
  charset c11_c31_n4l1 = , 1952-3490\3  4172-4951\3  9803-10090\3;
  charset c12_c22_c32_n22_n42_n4l2_n62 = , 1953-3490\3  3492-4171\3  4173-4951\3  7041-8080\3  8427-9802\3  9804-10090\3  11901-12421\3;
  charset c21 = , 3491-4171\3;
  charset cytb2_n12_n32_n41_n52 = , 4953-6085\3  6087-7039\3  8082-8425\3  8426-9802\3  10092-11899\3;
  charset n31 = , 8081-8425\3;
  charset n61 = , 11900-12421\3;
  charset n63 = , 11902-12421\3;
  charpartition mymodels =
    TN{17.6693,50.7353}+F{0.300639,0.331685,0.138115,0.229562}+G4{0.161794}: dloop{1.57476},
    HKY{28.7815}+F{0.396107,0.31873,0.1523,0.132863}+G4{0.221188}: a61_a81{0.139524},
    HKY{18.4013}+F{0.163742,0.305557,0.079646,0.451056}: a62{0.07087},
    TN{141.536,35.8405}+F{0.365027,0.432444,0.0491946,0.153334}+G4{0.899457}: a63_a83_c13_c23_c33_cytb3_n13_n23_n33_n43_n4l3_n53{0.218358},
    TN{64.3372,37.09}+F{0.376684,0.279998,0.145863,0.197455}+G4{0.0200244}: a82_cytb1_n11_n21_n51_12s_loops_16s_loops{0.0971929},
    K2P{292.724}+FQ+G4{0.0200235}: c11_c31_n4l1{0.0255592},
    HKY{24.3004}+F{0.192194,0.261047,0.132961,0.413798}+G4{0.0200385}: c12_c22_c32_n22_n42_n4l2_n62{0.0165576},
    K2P{399.88}+FQ+G4{0.020039}: c21{0.0690718},
    TPM3u{11.7886,56.1594}+F{0.232773,0.294339,0.118526,0.354362}+G4{0.0200288}: cytb2_n12_n32_n41_n52{0.039998},
    TNe{117.767,16.3359}+FQ+G4{0.0200391}: n31{0.216951},
    HKY{21.4725}+F{0.212098,0.0573749,0.466031,0.264496}+G4{0.0200237}: n61{0.0716608},
    TN{69.1058,199.979}+F{0.200989,0.0343543,0.350767,0.413889}: n63{0.231128};
end;

In principle (though it would take some checking, and converting file formats like this is notoriously tricky), one could borrow from the PartitionFinder code to write a translator to MrBayes. The relevant functions are:

• write_mrbayes(): https://github.com/brettc/partitionfinder/blob/42010508ea382c6450b114c32827ace8727a9d9d/partfinder/reporter.py#L236

• get_mrbayes_modeltext_dna(): https://github.com/brettc/partitionfinder/blob/42010508ea382c6450b114c32827ace8727a9d9d/partfinder/model_utils.py#L79

• get_mrbayes_modeltext_protein(): https://github.com/brettc/partitionfinder/blob/42010508ea382c6450b114c32827ace8727a9d9d/partfinder/model_utils.py#L110

Those functions will already do a lot of the heavy lifting in terms of figuring out what to convert to what, though it would be good to also look up if MrBayes has changed since ~2012 when I wrote them.

I don't think that this will make it close to the top of our list of things to do, but it would be worth adding a feature request to the issues, and referencing this discussion item. That way it will be there when we look, and if anyone in the community wants to write a script to do it, they can use this information as a jumping off point.

[bqminh]

We can request a new feature, that when -mset mrbayes option is specified, then IQ-TREE will generate a file ready to use for a subsequent MrBayes run. I don't know how many users would like this feature though. If there are at least 10, then Yes.

[roblanf]

@eosthusiast, if you want to go ahead and make a feature request (being as detailed as possible helps, especially if you can give a few examples, e.g. for DNA and for Protein models, suggest what should be done with models (like free rate models) which aren't in MrBayes, etc. etc.), and can then demonstrate enough interest as @bqminh suggests, then we will look into implementing this feature.

[protokute]

Oh, I think anyone who does partitioned analyses with MrBayes would benefit from this, as it will save some time, and some uncertainty from those not really familiar with the output files and the substitution model parameters.

[eosthusiast]

Thanks @bqminh, @roblanf and @protokute, I drew up the following feature request. Please let me know if there's anything I should add or change.

GitHub feature request: #195

[protokute]

Thanks a lot, @eosthusiast. I'd also add that it would be important that IQTREE also add all the specific relevant priors and parameters, the state frequency priors ("statefreqpr"), the reversible rate matrix prior ("revmatpr"), the gamma shape parameter ("shapepr") and the proportion of invariable sites parameter ("pinvar") when it is relevant for the selected substitution models.

This is an example of a MrBayes block that I have with all important priors and parameters:

begin mrbayes;

charset part1 = 1-269 427-691;
charset part2 = 270-426 1626-1807\3 2286-2356\3 3198-3521\3;
charset part3 = 692-1625 1628-1807\3 1870-2006\3 1871-2006\3 2072-2217\3 2288-2356\3 2898-3143\3 2899-3143\3 3197-3521\3;
charset part4 = 1627-1807\3 1872-2006\3 2071-2217\3 2287-2356\3 2897-3143\3;
charset part5 = 1808-1869 2007-2069;
charset part6 = 2070-2217\3 2357-2896;
charset part7 = 2218-2285 3144-3196 3199-3521\3;
partition Favolus = 7: part1, part2, part3, part4, part5, part6, part7;
set partition = Favolus;
outgroup Trametes_conchifer_FP106793sp;

lset applyto=(1) nst=2 rates=gamma ngammacat=4;
prset applyto=(1) revmatpr=fixed(1, 3.345048737, 1, 1, 3.345048737, 1) shapepr=fixed(0.477033) statefreqpr=dirichlet(0.224266,0.226241,0.237902,0.311591);

lset applyto=(2) nst=2 rates=propinv ngammacat=4;
prset applyto=(2) revmatpr=fixed(1, 9.142173815, 1, 1, 9.142173815, 1) pinvarpr=fixed(0.930324);

lset applyto=(3) nst=6 rates=invgamma ngammacat=4;
prset applyto=(3) revmatpr=fixed(1.299875016, 3.606095304, 1.996952633, 0.9918935694, 11.66345025, 1) statefreqpr=dirichlet(0.272392,0.205812,0.304472,0.217324) pinvarpr=fixed(0.700745);

lset applyto=(4) nst=2 rates=invgamma ngammacat=4;
prset applyto=(4) revmatpr=fixed(1, 7.965182054, 1, 1, 7.965182054, 1) shapepr=fixed(1.77888) statefreqpr=dirichlet(0.117995,0.396033,0.257864,0.228108) pinvarpr=fixed(0.195018);

lset applyto=(5) nst=2 rates=propinv ngammacat=4;
prset applyto=(5) revmatpr=fixed(1, 2.875739469, 1, 1, 2.875739469, 1) statefreqpr=dirichlet(0.222328,0.214478,0.206865,0.35633) pinvarpr=fixed(0.106448);

lset applyto=(6) nst=2 rates=gamma ngammacat=4;
prset applyto=(6) revmatpr=fixed(1, 4.363429921, 1, 1, 4.363429921, 1) statefreqpr=dirichlet(0.174006,0.290131,0.270597,0.265266) shapepr=fixed(0.368305);

lset applyto=(7) nst=6 rates=gamma ngammacat=4;
prset applyto=(7) revmatpr=fixed(4.672533274, 12.22736915, 2.36996032, 2.36696048, 13.49454344, 1) statefreqpr=dirichlet(0.182876,0.299937,0.221483,0.295704) shapepr=fixed(2.64497);

end;

[roblanf]

@protokute, this is useful.

Can you please add this to the feature request, but also give information about what it all means. I.e. what does fixed() imply, what do all the various commands refer to with respect to the output from IQ-TREE etc.

If we implement this, we need to understand precisely what people are looking for. None of us use MrBayes regularly (I don't think I've used it for >10 years).