generate and report partitions file now, update to v1.2.1

Author: AstrobioMike

bin/GToTree

@@ -4,7 +4,7 @@
 GREEN='\033[0;32m'
 RED='\033[0;31m'
 NC='\033[0m'
-VERSION="v1.1.13"
+VERSION="v1.2.1"
 
 
 if [ "$1" == "--version" ] || [ "$1" == "-v" ]; then
@@ -63,6 +63,9 @@ if [ "$#" == 0 ] || [ $1 == "-h" ] || [ $1 == "help" ]; then
     printf "                  the labels if using TaxonKit (-t flag specified). E.g., all would be\n"
     printf "                  \"-L Domain,Phylum,Class,Order,Family,Genus,Species,Strain\"\n\n"
 
+    printf "        - [-N] default: false\n"
+    printf "                  No tree. Generate alignment only.\n\n"    
+
     printf "        - [-T <str>] default: FastTree\n"
     printf "                  Which program to use for tree generation. Currently supported are\n"
     printf "                  \"FastTree\" and \"IQ-TREE\". As of now, these run with default settings\n"
@@ -149,9 +152,10 @@ taxonkit_id_swap='false'
 debug_flag='false'
 lineage_spec="Domain,Phylum,Class,Species,Strain"
 best_hit_mode='false'
+align_only='false'
 tree_program='FastTree'
 
-while getopts :a:g:f:A:H:o:m:tL:T:c:G:Bdn:j: args
+while getopts :a:g:f:A:H:o:m:tL:NT:c:G:Bdn:j: args
 do
     case "${args}"
     in
@@ -164,6 +168,7 @@ do
         m) file_to_genome_id_map=${OPTARG};;
         t) taxonkit_id_swap='true';;
         L) lineage_spec=${OPTARG};;
+        N) align_only='true' ;;
         T) tree_program=${OPTARG};;
         c) len_cutoff=${OPTARG};;
         G) gen_cutoff=${OPTARG};;
@@ -488,7 +493,7 @@ fi
 #############################################################################
 ##############  EXPLICITLY STATING IF DEFAULT BEHAVIOR CHANGED  #############
 #############################################################################
-if [ $output_dir != "GToTree_output" ] || [ "$file_to_genome_id_map" != "" ] || [ $taxonkit_id_swap != "false" ] || [ $len_cutoff != "0.2" ] || [ $gen_cutoff != "0.75" ] || [ $debug_flag == "true" ] || [ $best_hit_mode == "true" ] || [ $num_jobs != "1" ] || [ $num_cpus != 2 ] || [ $lineage_spec != "Domain,Phylum,Class,Species,Strain" ] || [ $tree_program != "FastTree" ]; then
+if [ $output_dir != "GToTree_output" ] || [ "$file_to_genome_id_map" != "" ] || [ $taxonkit_id_swap != "false" ] || [ $len_cutoff != "0.2" ] || [ $gen_cutoff != "0.75" ] || [ $debug_flag == "true" ] || [ $best_hit_mode == "true" ] || [ $num_jobs != "1" ] || [ $num_cpus != 2 ] || [ $lineage_spec != "Domain,Phylum,Class,Species,Strain" ] || [ $tree_program != "FastTree" ] || [ $align_only == 'true' ]; then
 
     if [ "$file_to_genome_id_map" != "" ]; then
         if [ ! -s $file_to_genome_id_map ]; then
@@ -521,8 +526,16 @@ if [ $output_dir != "GToTree_output" ] || [ "$file_to_genome_id_map" != "" ] ||
         printf "      - Gene-length filtering cutoff threshold (\"-c\") has been set to $len_cutoff.\n" | tee -a $gtotree_log
     fi
 
-    if [ $tree_program != "FastTree" ]; then
-        printf "      - Tree generation program (\"-T\") has been set to $tree_program.\n" | tee -a $gtotree_log
+    if [ $align_only == 'true' ]; then
+        printf "      - Only generating alignment, no tree, as \"-N\" option has been provided.\n" | tee -a $gtotree_log
+    fi
+
+    if [ $align_only == 'false' ]; then
+
+        if [ $tree_program != "FastTree" ]; then
+            printf "      - Tree generation program (\"-T\") has been set to $tree_program.\n" | tee -a $gtotree_log
+        fi
+
     fi
 
     if [ $gen_cutoff != "0.5" ]; then
@@ -1579,6 +1592,42 @@ grep ">" ${tmp_dir}/cat.tmp | cut -f1 > ${tmp_dir}/headers.tmp
 grep -v ">" ${tmp_dir}/cat.tmp | sed -e $'s/\t/XXXXX/g' > ${tmp_dir}/seqs.tmp
 paste -d "\n" ${tmp_dir}/headers.tmp ${tmp_dir}/seqs.tmp > ${output_dir}/Aligned_SCGs.faa
 
+# creating partitions file of gene coordinates for those that want might want it for mixed models with something like iqtree
+curr_start=1
+curr_stop=0
+n=1
+num_to_do=$(ls ${tmp_dir}/*_all_aligned.faa | wc -l)
+
+for gene_file in $(ls ${tmp_dir}/*_all_aligned.faa)
+
+do
+    if [ $n == 1 ]; then
+        n=$(($n + 1))
+        gene_base=$(basename $gene_file)
+        gene=${gene_base%%_all_aligned.faa}
+        len_tmp=$(head -n 2 $gene_file | grep -v ">" | wc -c)
+        len_tmp=$(echo "$len_tmp - 1" | bc)
+
+        printf "AA, $gene = 1-${len_tmp}\n"
+
+        curr_start=$(echo "$len_tmp + 6" | bc)
+
+    else
+        n=$(($n + 1))
+        gene_base=$(basename $gene_file)
+        gene=${gene_base%%_all_aligned.faa}
+        len_tmp=$(head -n 2 $gene_file | grep -v ">" | wc -c)
+        len_tmp=$(echo "$len_tmp - 1" | bc)
+
+        curr_stop=$(echo "$curr_start + $len_tmp - 1" | bc)
+
+        printf "AA, $gene = ${curr_start}-${curr_stop}\n"
+
+        curr_start=$(echo "$curr_stop + 6" | bc)
+    fi
+
+done > ${output_dir}/Partitions.txt
+
 # storing genomes that made it through workflow to report at end
 genomes_retained=$(wc -l ${tmp_dir}/final_genomes_from_all_sources.tmp | sed "s/^ *//" | cut -d " " -f 1)
 
@@ -1938,47 +1987,51 @@ fi
 ##############################   MAKING TREE  ###############################
 #############################################################################
 
-if [ $tree_program == 'FastTree' ]; then
+if [ $align_only == 'false' ]; then
 
-    printf "\n ############################################################################## \n" | tee -a $gtotree_log
-    printf " ####                          Running FastTree                            ####\n" | tee -a $gtotree_log
-    printf " ############################################################################## \n\n" | tee -a $gtotree_log
+    if [ $tree_program == 'FastTree' ]; then
 
-    curr_time=$(date +"%I:%M %p")
-    duration=$SECONDS
+        printf "\n ############################################################################## \n" | tee -a $gtotree_log
+        printf " ####                          Running FastTree                            ####\n" | tee -a $gtotree_log
+        printf " ############################################################################## \n\n" | tee -a $gtotree_log
 
-    printf "           It is currently $curr_time; the process started at $start_time.\n" | tee -a $gtotree_log
-    printf "               Current process runtime: $(($duration / 60 / 60)) hours and $((($duration / 60) % 60)) minutes.\n\n" | tee -a $gtotree_log
+        curr_time=$(date +"%I:%M %p")
+        duration=$SECONDS
 
+        printf "           It is currently $curr_time; the process started at $start_time.\n" | tee -a $gtotree_log
+        printf "               Current process runtime: $(($duration / 60 / 60)) hours and $((($duration / 60) % 60)) minutes.\n\n" | tee -a $gtotree_log
+
+
+        if [ -s ${output_dir}/Aligned_SCGs_mod_names.faa ]; then
+            FastTree ${output_dir}/Aligned_SCGs_mod_names.faa > ${output_dir}/Aligned_SCGs_mod_names.tre | tee -a $gtotree_log
+        else
+            FastTree ${output_dir}/Aligned_SCGs.faa > ${output_dir}/Aligned_SCGs.tre | tee -a $gtotree_log
+        fi
 
-    if [ -s ${output_dir}/Aligned_SCGs_mod_names.faa ]; then
-        FastTree ${output_dir}/Aligned_SCGs_mod_names.faa > ${output_dir}/Aligned_SCGs_mod_names.tre | tee -a $gtotree_log
     else
-        FastTree ${output_dir}/Aligned_SCGs.faa > ${output_dir}/Aligned_SCGs.tre | tee -a $gtotree_log
-    fi
 
-else
+        printf "\n ############################################################################## \n" | tee -a $gtotree_log
+        printf " ####                          Running IQ-TREE                            ####\n" | tee -a $gtotree_log
+        printf " ############################################################################## \n\n" | tee -a $gtotree_log
 
-    printf "\n ############################################################################## \n" | tee -a $gtotree_log
-    printf " ####                          Running IQ-TREE                            ####\n" | tee -a $gtotree_log
-    printf " ############################################################################## \n\n" | tee -a $gtotree_log
+        curr_time=$(date +"%I:%M %p")
+        duration=$SECONDS
 
-    curr_time=$(date +"%I:%M %p")
-    duration=$SECONDS
+        printf "           It is currently $curr_time; the process started at $start_time.\n" | tee -a $gtotree_log
+        printf "               Current process runtime: $(($duration / 60 / 60)) hours and $((($duration / 60) % 60)) minutes.\n\n" | tee -a $gtotree_log
 
-    printf "           It is currently $curr_time; the process started at $start_time.\n" | tee -a $gtotree_log
-    printf "               Current process runtime: $(($duration / 60 / 60)) hours and $((($duration / 60) % 60)) minutes.\n\n" | tee -a $gtotree_log
 
+        if [ -s ${output_dir}/Aligned_SCGs_mod_names.faa ]; then
+            iqtree -s ${output_dir}/Aligned_SCGs_mod_names.faa -nt $num_jobs -mset WAG,LG -bb 1000 -pre iqtree_out
+        else
+            iqtree -s ${output_dir}/Aligned_SCGs.faa -nt $num_jobs -mset WAG,LG -bb 1000 -pre iqtree_out
+        fi
 
-    if [ -s ${output_dir}/Aligned_SCGs_mod_names.faa ]; then
-        iqtree -s ${output_dir}/Aligned_SCGs_mod_names.faa -nt $num_jobs -mset WAG,LG -bb 1000 -pre iqtree_out
-    else
-        iqtree -s ${output_dir}/Aligned_SCGs.faa -nt $num_jobs -mset WAG,LG -bb 1000 -pre iqtree_out
-    fi
+        mkdir ${output_dir}/iqtree_out/
+        mv iqtree_out* ${output_dir}/iqtree_out/
+        cp ${output_dir}/iqtree_out/iqtree_out.treefile ${output_dir}/Aligned_SCGs_mod_names.tre
 
-    mkdir ${output_dir}/iqtree_out/
-    mv iqtree_out* ${output_dir}/iqtree_out/
-    cp ${output_dir}/iqtree_out/iqtree_out.treefile ${output_dir}/Aligned_SCGs_mod_names.tre
+    fi
 
 fi
 
@@ -2000,22 +2053,34 @@ fi
 
 # reporting primary output files
 
-printf "    Full alignment written to file:\n" | tee -a $gtotree_log
+if [ $align_only == 'false' ]; then
+
+    if [ -s ${output_dir}/Aligned_SCGs_mod_names.faa ]; then
+
+        printf "    Tree written to:\n" | tee -a $gtotree_log
+        printf "        ${GREEN}${output_dir}/Aligned_SCGs_mod_names.tre${NC}\n\n" | tee -a $gtotree_log
+    else
+        printf "    Tree written to:\n" | tee -a $gtotree_log
+        printf "        ${GREEN}${output_dir}/Aligned_SCGs.tre${NC}\n\n" | tee -a $gtotree_log
+    fi
+
+fi
+
+printf "    Full alignment written to:\n" | tee -a $gtotree_log
 printf "        ${GREEN}${output_dir}/Aligned_SCGs.faa${NC}\n\n" | tee -a $gtotree_log
+
 if [ -s ${output_dir}/Aligned_SCGs_mod_names.faa ]; then
-    printf "    Alignment with altered headers written to file:\n" | tee -a $gtotree_log
+    printf "    Alignment with altered headers written to:\n" | tee -a $gtotree_log
     printf "        ${GREEN}${output_dir}/Aligned_SCGs_mod_names.faa${NC}\n\n" | tee -a $gtotree_log
-    printf "    Tree written to file:\n" | tee -a $gtotree_log
-    printf "        ${GREEN}${output_dir}/Aligned_SCGs_mod_names.tre${NC}\n\n" | tee -a $gtotree_log
-else
-    printf "    Tree written to file:\n" | tee -a $gtotree_log
-    printf "        ${GREEN}${output_dir}/Aligned_SCGs.tre${NC}\n\n" | tee -a $gtotree_log
 fi
 
-printf "    Summary file with comp./redund. estimates written to file:\n" | tee -a $gtotree_log
+printf "    Partitions (for downstream use with mixed-model treeing) written to:\n" | tee -a $gtotree_log
+printf "        ${GREEN}${output_dir}/Partitions.txt${NC}\n\n" | tee -a $gtotree_log
+
+printf "    Summary file with comp./redund. estimates written to:\n" | tee -a $gtotree_log
 printf "        ${GREEN}${output_dir}/All_genomes_summary_info.tsv${NC}\n\n" | tee -a $gtotree_log
 
-printf "    Summary table with hits per target gene per genome written to file:\n" | tee -a $gtotree_log
+printf "    Summary table with hits per target gene per genome written to:\n" | tee -a $gtotree_log
 printf "        ${GREEN}${output_dir}/All_genomes_SCG_hit_counts.tsv${NC}\n\n" | tee -a $gtotree_log
 
 # reporting any problem files/accessions

example_run/Alteromonas_example/Aligned_SCGs.faa

@@ -1 +1 @@
-((GCF_000020585.3_Alteromonas_mediterranea_DE:0.00202,(GCF_000439535.1_Alteromonas_mediterranea_U7:0.0,GCF_000439555.1_Alteromonas_mediterranea_U8:0.0):0.00153)0.900:0.00055,((GCF_001886395.1_Alteromonas_mediterranea_CP48:0.00138,(GCF_001886435.1_Alteromonas_mediterranea_RG65:0.00132,(GCF_001886415.1_Alteromonas_mediterranea_CP49:0.00055,GCF_001886455.1_Alteromonas_mediterranea_AR43:0.00055)1.000:0.00189)0.881:0.00050)1.000:0.00083,(((GCF_001886475.1_Alteromonas_sp_RW2A1_RW2A1:0.00124,GCF_002831605.1_Alteromonas_sp_MB-3u-76_MB-3u-76:0.00094)1.000:0.08194,((GCF_000730385.1_Alteromonas_australica_H_17:0.00122,GCF_000934525.1_Alteromonas_australica_DE170:0.00138)1.000:0.06762,(GCF_000011365.1_Alpha_Outgroup:1.15277,(GCF_000213655.1_Alteromonas_naphthalenivorans_SN2:0.00831,((GCF_001433715.1_Alteromonas_stellipolaris_LMG_21856:0.0,GCF_001562115.1_Alteromonas_stellipolaris_LMG_21861:0.0):0.00159,(GCF_001632825.1_Alteromonas_stellipolaris_PQQ-44:0.00080,(GCF_001632765.1_Alteromonas_stellipolaris_PQQ-42:0.00055,(GCF_001562125.1_Alteromonas_sp_Mac1_Mac1:0.00114,GCF_001562195.1_Alteromonas_addita_R10SW13:0.00109)0.872:0.00051)0.347:0.00055)0.999:0.00104)1.000:0.00810)1.000:0.05790)0.982:0.01987)1.000:0.03197)1.000:0.03452,(GCF_003443615.1_Alteromonas_sp_BL110_BL110:0.01952,(GCF_000299995.1_Alteromonas_macleodii_Black_Sea_11:0.01052,(((((Our_Alteromonas_MAG:0.00221,GCF_002849875.1_Alteromonas_macleodii_Te101:0.00182)0.643:0.00055,(GCF_000172635.2_Alteromonas_macleodii_ATCC_27126:0.00146,GCF_000299955.1_Alteromonas_macleodii_English_Channel_673:0.00215)0.882:0.00055)0.670:0.00055,GCF_001578515.1_Alteromonas_macleodii_HOT1A3:0.00198)0.999:0.00059,GCF_000300175.1_Alteromonas_macleodii_Balearic_Sea_AD45:0.00364)0.867:0.00052,GCF_001562235.1_Alteromonas_macleodii_D7:0.00364)1.000:0.00561)1.000:0.01088)1.000:0.01130)1.000:0.02062)1.000:0.00084,(GCF_000439495.1_Alteromonas_mediterranea_MED64:0.00126,((GCF_000310085.1_Alteromonas_mediterranea_AltDE1:0.0,GCF_000439575.1_Alteromonas_mediterranea_UM7:0.0,GCF_001562315.1_Alteromonas_mediterranea_UM8:0.0):0.00055,GCF_000439595.1_Alteromonas_mediterranea_UM4b:0.00055)1.000:0.00166)0.059:0.00055);
+((GCF_000020585.3_Alteromonas_mediterranea_DE:0.00207,GCF_000439495.1_Alteromonas_mediterranea_MED64:0.00119)0.724:0.00055,((GCF_001886395.1_Alteromonas_mediterranea_CP48:0.00135,(GCF_001886435.1_Alteromonas_mediterranea_RG65:0.00126,(GCF_001886415.1_Alteromonas_mediterranea_CP49:0.00055,GCF_001886455.1_Alteromonas_mediterranea_AR43:0.00055)1.000:0.00185)0.916:0.00055)1.000:0.00081,(((GCF_001886475.1_Alteromonas_sp_RW2A1_RW2A1:0.00124,GCF_002831605.1_Alteromonas_sp_MB-3u-76_MB-3u-76:0.00094)1.000:0.08171,((GCF_000730385.1_Alteromonas_australica_H_17:0.00118,GCF_000934525.1_Alteromonas_australica_DE170:0.00140)1.000:0.06761,(GCF_000011365.1_Alpha_Outgroup:1.15183,(GCF_000213655.1_Alteromonas_naphthalenivorans_SN2:0.00828,((GCF_001433715.1_Alteromonas_stellipolaris_LMG_21856:0.0,GCF_001562115.1_Alteromonas_stellipolaris_LMG_21861:0.0):0.00150,(GCF_001632825.1_Alteromonas_stellipolaris_PQQ-44:0.00080,(GCF_001632765.1_Alteromonas_stellipolaris_PQQ-42:0.00055,(GCF_001562125.1_Alteromonas_sp_Mac1_Mac1:0.00114,GCF_001562195.1_Alteromonas_addita_R10SW13:0.00110)0.785:0.00050)0.317:0.00055)0.998:0.00104)1.000:0.00824)1.000:0.05757)0.986:0.02015)1.000:0.03187)1.000:0.03449,(GCF_003443615.1_Alteromonas_sp_BL110_BL110:0.01934,(GCF_000299995.1_Alteromonas_macleodii_Black_Sea_11:0.01022,(GCF_001562235.1_Alteromonas_macleodii_D7:0.00361,(GCF_000300175.1_Alteromonas_macleodii_Balearic_Sea_AD45:0.00363,(GCF_001578515.1_Alteromonas_macleodii_HOT1A3:0.00198,((GCF_000172635.2_Alteromonas_macleodii_ATCC_27126:0.00146,GCF_000299955.1_Alteromonas_macleodii_English_Channel_673:0.00215)0.875:0.00055,(GCF_002849875.1_Alteromonas_macleodii_Te101:0.00182,Our_Alteromonas_MAG:0.00214)0.601:0.00055)0.657:0.00055)1.000:0.00059)0.863:0.00050)1.000:0.00535)1.000:0.01076)1.000:0.01117)1.000:0.02013)1.000:0.00077,(((GCF_000310085.1_Alteromonas_mediterranea_AltDE1:0.0,GCF_000439575.1_Alteromonas_mediterranea_UM7:0.0,GCF_001562315.1_Alteromonas_mediterranea_UM8:0.0):0.00055,GCF_000439595.1_Alteromonas_mediterranea_UM4b:0.00055)1.000:0.00152,(GCF_000439535.1_Alteromonas_mediterranea_U7:0.0,GCF_000439555.1_Alteromonas_mediterranea_U8:0.0):0.00167)0.924:0.00055);