From 37c1c36c7ab14ff0441e13caae6e5f5aacfef6d0 Mon Sep 17 00:00:00 2001 From: Auden Cote-L'Heureux <52716489+AudenCote@users.noreply.github.com> Date: Mon, 21 Aug 2023 13:07:07 -0400 Subject: [PATCH] Adding date to stat summary script --- PTL1/Genomes/Scripts/6_SummaryStats.py | 27 ++++++++++++++------------ 1 file changed, 15 insertions(+), 12 deletions(-) diff --git a/PTL1/Genomes/Scripts/6_SummaryStats.py b/PTL1/Genomes/Scripts/6_SummaryStats.py index f892e2e..e3d4967 100644 --- a/PTL1/Genomes/Scripts/6_SummaryStats.py +++ b/PTL1/Genomes/Scripts/6_SummaryStats.py @@ -7,8 +7,11 @@ from math import ceil, floor from tqdm import tqdm #import matplotlib.pyplot as plt import numpy as np +from datetime import date +today = date.today() + def get_args(): parser = argparse.ArgumentParser( @@ -102,13 +105,13 @@ def get_nuc_comp(args, gcodes): def per_seq(args, nuc_comp, aa_comp, all_transcripts, r2g_lengths, transcript_id_corr, og_mean_lens): - if not os.path.isdir(args.input + '/PerSequenceStatSummaries'): - os.mkdir(args.input + '/PerSequenceStatSummaries') + if not os.path.isdir(args.input + '/PerSequenceStatSummaries_' + today): + os.mkdir(args.input + '/PerSequenceStatSummaries_' + today) taxa = list(dict.fromkeys([seq[:10] for seq in nuc_comp])) for taxon in taxa: - with open(args.input + '/PerSequenceStatSummaries/' + taxon + '.csv', 'w') as o: + with open(args.input + '/PerSequenceStatSummaries_' + today + '/' + taxon + '.csv', 'w') as o: o.write('Sequence,Taxon,OG,OrigName,OrigLength,R2GLength,AvgLengthOGinHook,AmbiguousCodons,GC-Overall,GC1,GC2,GC3,GC3-Degen,ExpWrightENc,ObsWrightENc_6Fold,ObsWrightENc_No6Fold,ObsWeightedENc_6Fold,ObsWeightedENc_No6Fold,FYMINK,GARP,OtherAA,N.Xs\n') for rec in nuc_comp: if rec[:10] == taxon: @@ -133,7 +136,7 @@ def per_tax(args, nuc_comp, aa_comp, all_transcripts, r2g_lengths, gcodes, og_me taxa = list(dict.fromkeys([seq[:10] for seq in nuc_comp])) - with open(args.input + '/PerTaxonSummary.csv', 'w') as o: + with open(args.input + '/PerTaxonSummary_' + today + '.csv', 'w') as o: o.write('Taxon,OrigSeqs,Orig_MedianGC,Orig_GCWidth_5-95Perc,Orig_MedianLen,Orig_IQRLen,R2GSeqs,R2GOGs,R2GMedian_GC3,R2G_5Perc_GC3,R2G_95Perc_GC3,R2G_GC3Width_5-95Perc,R2G_MedianENc,R2G_IQRENc,R2G_MedianLen,R2G_IQRLen,R2G_Prop.G1.5_OGAvg,R2G_Prop.L0.5_OGAvg,R2G_MeanXs,GeneticCode\n') for taxon in taxa: @@ -190,11 +193,11 @@ def r2g_jf(args, nuc_comp, gcodes): #Q: should there be an maximum IQR cutoff at which we do NOT produce a file here? - if not os.path.isdir(args.input + '/ReadyToGo/ReadyToGo_NTD_JF'): - os.mkdir(args.input + '/ReadyToGo/ReadyToGo_NTD_JF') + if not os.path.isdir(args.input + '/ReadyToGo/ReadyToGo_NTD_JF_' + today): + os.mkdir(args.input + '/ReadyToGo/ReadyToGo_NTD_JF_' + today) - if not os.path.isdir(args.input + '/ReadyToGo/ReadyToGo_AA_JF'): - os.mkdir(args.input + '/ReadyToGo/ReadyToGo_AA_JF') + if not os.path.isdir(args.input + '/ReadyToGo/ReadyToGo_AA_JF_' + today): + os.mkdir(args.input + '/ReadyToGo/ReadyToGo_AA_JF_' + today) for file in os.listdir(args.input + '/ReadyToGo/ReadyToGo_NTD'): if file.endswith('.fasta') and file[:10] in gcodes: @@ -203,12 +206,12 @@ def r2g_jf(args, nuc_comp, gcodes): r2g_ntds = [nuc_comp[seq] for seq in nuc_comp if seq[:10] == taxon] r2g_gc3s = sorted([seq.gc4F for seq in r2g_ntds]) - with open(args.input + '/ReadyToGo/ReadyToGo_NTD_JF/' + file.replace('.fasta', '.JF.fasta'), 'w') as o: + with open(args.input + '/ReadyToGo/ReadyToGo_NTD_JF_' + today + '/' + file.replace('.fasta', '.JF.fasta'), 'w') as o: for rec in SeqIO.parse(args.input + '/ReadyToGo/ReadyToGo_NTD/' + file, 'fasta'): if nuc_comp[rec.id].gc4F > r2g_gc3s[floor(len(r2g_gc3s)*0.25)] and nuc_comp[rec.id].gc4F < r2g_gc3s[floor(len(r2g_gc3s)*0.75)]: o.write('>' + rec.id + '\n' + str(rec.seq) + '\n\n') - with open(args.input + '/ReadyToGo/ReadyToGo_AA_JF/' + file.replace('.fasta', '.JF.fasta').replace('NTD', 'AA'), 'w') as o: + with open(args.input + '/ReadyToGo/ReadyToGo_AA_JF_' + today + '/' + file.replace('.fasta', '.JF.fasta').replace('NTD', 'AA'), 'w') as o: for rec in SeqIO.parse(args.input + '/ReadyToGo/ReadyToGo_AA/' + file.replace('NTD', 'AA'), 'fasta'): if nuc_comp[rec.id].gc4F > r2g_gc3s[floor(len(r2g_gc3s)*0.25)] and nuc_comp[rec.id].gc4F < r2g_gc3s[floor(len(r2g_gc3s)*0.75)]: o.write('>' + rec.id + '\n' + str(rec.seq) + '\n\n') @@ -216,8 +219,8 @@ def r2g_jf(args, nuc_comp, gcodes): def plot_jf(args, nuc_comp): - if not os.path.isdir(args.input + '/GC3xENc_Plots'): - os.mkdir(args.input + '/GC3xENc_Plots') + if not os.path.isdir(args.input + '/GC3xENc_Plots_' + today): + os.mkdir(args.input + '/GC3xENc_Plots_' + today) taxa = list(dict.fromkeys([rec[:10] for rec in nuc_comp]))