#Code for Quattrini et al. Mito-nuclear discordance in Anthozoa, with notes on unique properties of their mitochondrial genomes #Follow Phyluce tutorial #https://phyluce.readthedocs.io/en/latest/tutorials/tutorial-1.html with modifications described in Quattrini et al. 20218, 2020 phyluce_assembly_get_match_counts --locus-db uce-search-results/probe.matches.sqlite --taxon-list-config taxonset.conf --taxon-group 'all' --incomplete-matrix --output taxon-sets/all/all-uce-taxa-incomplete.conf phyluce_assembly_get_match_counts --locus-db trans-search-results/probe.matches.sqlite --taxon-list-config taxonset.conf --taxon-group 'all' --incomplete-matrix --output taxon-sets/all/all-trans-taxa-incomplete.conf phyluce_assembly_get_fastas_from_match_counts --contigs /scratch/nmnh_corals/quattrinia/analysis/UCE/mito_nucl/phyluce-octos/contigs --locus-db /scratch/nmnh_corals/quattrinia/analysis/UCE/mito_nucl/phyluce-octos/trans-search-results/probe.matches.sqlite --match-count-output all-trans-taxa-incomplete.conf --output all-trans-incomplete.fasta --incomplete-matrix all-trans-taxa-incomplete.incomplete --log-path log phyluce_assembly_get_fastas_from_match_counts --contigs /scratch/nmnh_corals/quattrinia/analysis/UCE/mito_nucl/phyluce-octos/contigs --locus-db /scratch/nmnh_corals/quattrinia/analysis/UCE/mito_nucl/phyluce-octos/uce-search-results/probe.matches.sqlite --match-count-output all-uce-taxa-incomplete.conf --output all-uce-incomplete.fasta --incomplete-matrix all-uce-taxa-incomplete.incomplete --log-path log phyluce_align_seqcap_align --input /scratch/nmnh_corals/quattrinia/analysis/UCE/mito_nucl/phyluce/taxon-sets/all/all-combin-trans.final.fasta --output /scratch/nmnh_corals/quattrinia/analysis/UCE/mito_nucl/phyluce/taxon-sets/all/mafft-trans-nexus-edge-trimmed --taxa 108 --aligner mafft --cores $NSLOTS --incomplete-matrix --log-path /scratch/nmnh_corals/quattrinia/analysis/UCE/mito_nucl/phyluce/taxon-sets/all/log phyluce_align_seqcap_align --input /scratch/nmnh_corals/quattrinia/analysis/UCE/mito_nucl/phyluce/taxon-sets/all/all-combin-uce.final.fasta --output /scratch/nmnh_corals/quattrinia/analysis/UCE/mito_nucl/phyluce/taxon-sets/all/mafft-uce-nexus-edge-trimmed --taxa 108 --aligner mafft --cores $NSLOTS --incomplete-matrix --log-path /scratch/nmnh_corals/quattrinia/analysis/UCE/mito_nucl/phyluce/taxon-sets/all/log phyluce_align_seqcap_align --input /scratch/nmnh_corals/quattrinia/analysis/UCE/mito_nucl/phyluce/taxon-sets/all/all-combin-trans.final.fasta --output /scratch/nmnh_corals/quattrinia/analysis/UCE/mito_nucl/phyluce/taxon-sets/all/mafft-trans-nexus-edge-trimmed --taxa 108 --aligner mafft --cores $NSLOTS --incomplete-matrix --log-path /scratch/nmnh_corals/quattrinia/analysis/UCE/mito_nucl/phyluce/taxon-sets/all/log --no-trim --output-format fasta phyluce_align_get_gblocks_trimmed_alignments_from_untrimmed --alignments /scratch/nmnh_corals/quattrinia/analysis/UCE/mito_nucl/phyluce/taxon-sets/all/mafft-uce-nexus-edge-trimmed --output /scratch/nmnh_corals/quattrinia/analysis/UCE/mito_nucl/phyluce/taxon-sets/all/mafft-uce-nexus-edge-trimmed-gblocks --cores $NSLOTS --log log --b1 0.5 --b2 0.5 --b3 10 --b4 5 phyluce_align_get_gblocks_trimmed_alignments_from_untrimmed --alignments /scratch/nmnh_corals/quattrinia/analysis/UCE/mito_nucl/phyluce/taxon-sets/all/mafft-trans-nexus-edge-trimmed --output /scratch/nmnh_corals/quattrinia/analysis/UCE/mito_nucl/phyluce/taxon-sets/all/mafft-trans-nexus-edge-trimmed-gblocks --cores $NSLOTS --log log --b1 0.5 --b2 0.5 --b3 10 --b4 5 phyluce_align_get_align_summary_data --alignments mafft-uce-nexus-edge-trimmed-gblocks --cores 12 --log-path log phyluce_align_get_align_summary_data --alignments mafft-trans-nexus-edge-trimmed-gblocks --cores 12 --log-path log #cp trans and uce loci together in one directory phyluce_align_remove_locus_name_from_files --alignments mafft-all-nexus-internal-trimmed-gblocks --output mafft-all-nexus-internal-trimmed-gblocks-clean --cores 12 --log-path log phyluce_align_get_only_loci_with_min_taxa --alignments mafft-all-nexus-internal-trimmed-gblocks-clean --taxa 94 --percent 0.75 --output mafft-all-nexus-internal-trimmed-gblocks-clean-75p --cores 12 --log-path log phyluce_align_concatenate_alignments --alignments mafft-all-nexus-internal-trimmed-gblocks-clean-75p --output mafft-all-nexus-internal-trimmed-gblocks-clean-75p-raxml --phylip --log-path log phyluce_align_concatenate_alignments --alignments mafft-all-nexus-internal-trimmed-gblocks-clean-50p --output mafft-all-nexus-internal-trimmed-gblocks-clean-50p-raxml --phylip --log-path log #run iqtree2 example iqtree2 -s /scratch/nmnh_corals/quattrinia/analysis/UCE/mito_nucl/hexa_nooutg_march22-concat/hexa_nooutg_march22-concat.phylip -p /scratch/nmnh_corals/quattrinia/analysis/UCE/mito_nucl/hexa_nooutg_march22-concat/hexa_nooutg_march22-concat.charsets -m TESTMERGE -pre hexa-nooutg-mar22 -rcluster 10 -bb 1000 -T $NSLOTS -undo -alrt 1000 ##R stuff #Drop R. muelleri from tree library(ape) library(ggtree) mito=read.tree("/Users/quattrinia/Manuscripts_and_Reports/UCES_mtgenomes/finaltrees/octo_75p_mar22.withlabs.tre") species=c("ANT111_Renilla_muelleri") mito<-drop.tip(mito,mito$tip.label[match(species, mito$tip.label)]) t50=read.tree("/Users/quattrinia/Manuscripts_and_Reports/UCES_mtgenomes/finaltrees/octo_50p_mar22.withlabs.tre") write.tree(mito) #to compare tip names etc between two trees comparePhylo(t50,mitopb) #treeshr and astral #run treeshrink to remove long branches rename .treefile .tree *.treefile cat *tree >out.trees /scratch/nmnh_corals/programs/TreeShrink/run_treeshrink.py -t out.trees -o treeshr java -jar /scratch/nmnh_corals/programs/ASTRAL/astral.5.7.8.jar -i output.trees -o astral.75lr #GC-content seqkit fx2tab --name --only-id --gc ANT*.fasta #average across loci for each sample awk -v N=2 '{ sum += $N } END { if (NR > 0) print sum / NR }' ANT15.gc.txt >ANT15.gc.avg.txt #Robinson-Foulds- Was run in the option of iqtree #concatenate all trees together in one file iqtree2 -rf_all octo.a.50lr.50lrm.50.75lr.75lrm.75.mt.mtpb.trees #TreeDist library(adegenet) library(TreeDist) library(ggtree) mytrees=read.tree("/Users/quattrinia/Manuscripts_and_Reports/UCES_mtgenomes/hexa.a.50lr.50lrm.50.75lr.75lrm.75.mt.mtpb.trees") TreeDistance(mytrees) #BoxPlots in R library(ggplot2) gctable=read.table("/Users/quattrinia/Manuscripts_and_Reports/UCES_mtgenomes/gc_content.table.txt") p=ggplot(gctable, aes(x=V1, y=V3, fill=V2)) + + geom_boxplot()+ ylab("GC Content") + xlab("type") plot(p) #anova in R > gcocttable=read.table("/Users/quattrinia/Manuscripts_and_Reports/UCES_mtgenomes/octo_gc.txt") > gcoctaov=aov(V3~V2, data=gcocttable) Df Sum Sq Mean Sq F value Pr(>F) V2 1 1623.6 1623.6 778.6 <2e-16 *** Residuals 186 387.9 2.1 --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 > gchextable=read.table("/Users/quattrinia/Manuscripts_and_Reports/UCES_mtgenomes/hexa_gc.txt") > gchexaov=aov(V3~V1, data=gchextable) > summary(gchexaov) ummary(gchexaov) Df Sum Sq Mean Sq F value Pr(>F) V2 1 1835 1835.4 246.2 <2e-16 *** Residuals 214 1596 7.5 --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1