added output parameters

Author: aim11

README.md

@@ -60,31 +60,29 @@ Options:
 
 ### Examples
 Local Annotator
-Annotate somatic variants and CNAs using the local annotator:
+1. Annotate somatic variants and CNAs using the local annotator:
 > python main.py --annotator local --output path/to/output --somatic_variants path/to/snvs.tsv --ascatestimates path/to/ascat.tsv
 
-Annotate CNAs using the local annotator:
+2. Annotate CNAs using the local annotator:
 > python main.py --annotator local --output path/to/output --copy_number_alterations path/to/cnas.tsv --ascatestimates path/to/ascat.tsv
 
-External Annotator
-Annotate CNAs using OncoKB:
-> python main.py --annotator external --output path/to/output --oncokbcna --copy_number_alterations path/to/cnas.tsv
+External Annotator (execution order is important as the CGI annotations are supplementing the OncoKB annotations in default)
+3. Annotate CNAs using OncoKB:
+> python main.py --annotator external --output path/to/output --oncokbcna --copy_number_alterations path/to/locally_annotated_cnas.tsv
 
-Annotate somatic variants using OncoKB:
->python main.py --annotator external --output path/to/output --oncokbsnv --somatic_variants path/to/snvs.tsv
+4. Annotate somatic variants using OncoKB:
+>python main.py --annotator external --output path/to/output --oncokbsnv --somatic_variants path/to/locally_annotated_snvs.tsv
 
-Annotate CNAs using Cancer Genome Interpreter:
+5. Annotate CNAs using Cancer Genome Interpreter:
 
->python external_annotator.py --cgiquery --copy_number_alterations path/to/cnas.tsv
+>python external_annotator.py --cgiquery --copy_number_alterations path/to/oncokb_annotated_cnas.tsv
 
-Annotate somatic variants using Cancer Genome Interpreter:
->python main.py --annotator external --output path/to/output --cgiquery --somatic_variants path/to/snvs.tsv
+6. Annotate somatic variants using Cancer Genome Interpreter:
+>python main.py --annotator external --output path/to/output --cgiquery --somatic_variants path/to/oncokb_annotated_snvs.tsv
 
-SLURM Scripts
+SLURM Scripts: edit the scripts to set the correct paths and SLURM sbatch parameters.
 
 Submit a batch job to SLURM cluster to annotate on multiple computing nodes:
 >./slurm_scripts/annotate_cnas.sh path/to/sample_list.txt
 
 >./slurm_scripts/snv_annotation.sbatch path/to/sample_list.txt
-
-### License

cgi_annotator.py

@@ -105,7 +105,7 @@ def generate_cgi_cna_file_from_list(genelist):
             file2.write(row)
         file2.close()
 
-def launch_cgi_job_with_mulitple_variant_types(mutations_file, cnas_file, transloc_file, cancer_type, reference):
+def launch_cgi_job_with_mulitple_variant_types(mutations_file=None, cnas_file=None, transloc_file=None, cancer_type="HGSOC", reference="GRCh38"):
     """
         This function launches a CGI (Cancer Genome Interpreter) job with multiple variant types,
         using the CGI API. It takes in mutation, cnas, and translocation files, cancer type, and
@@ -174,7 +174,7 @@ def launch_cgi_job_with_mulitple_variant_types(mutations_file, cnas_file, transl
         return 0
 
 
-def query_cgi_job(jobid, snv_annotations: pd.DataFrame = None, cna_annotations: pd.DataFrame = None):
+def query_cgi_job(jobid, output, snv_annotations: pd.DataFrame = None, cna_annotations: pd.DataFrame = None):
     """
     Query the CGI API with a job ID and save the results to the database.
 
@@ -206,10 +206,8 @@ def query_cgi_job(jobid, snv_annotations: pd.DataFrame = None, cna_annotations:
         cgi_snvdf = None
         cgi_cnadf = None
         treatments = []
+
         for fn in fnames:
-            # reader = z.open(f)
-            # for row in reader.readlines():
-            #    print(row)
             z.extract(fn)
             df = pd.read_csv(fn, sep="\t")
             print(fn)
@@ -277,14 +275,14 @@ def query_cgi_job(jobid, snv_annotations: pd.DataFrame = None, cna_annotations:
                     snv_annotations.at[indxs, 'tumorTypeSummary'] = handle_string_field(cgi_snv["driver_statement"])
 
         if isinstance(snv_annotations, pd.DataFrame):
-            snv_annotations.to_csv("snv_annotated_cgi.csv", index=False, sep="\t", columns=['patient_id', 'sample_id', 'alteration', 'hugoSymbol', 'tumorType', 'consequence', 'oncogenic', 'mutationEffectDescription', 'gene_role', 'citationPMids', 'level_of_evidence', 'geneSummary', 'variantSummary', 'tumorTypeSummary'])
+            snv_annotations.to_csv(output, index=False, sep="\t", columns=['patient_id', 'sample_id', 'alteration', 'hugoSymbol', 'tumorType', 'consequence', 'oncogenic', 'mutationEffectDescription', 'gene_role', 'citationPMids', 'level_of_evidence', 'geneSummary', 'variantSummary', 'tumorTypeSummary'])
             trdf = pd.DataFrame(treatments)
-            trdf.to_csv("treatments_cgi_snv.csv", index=False, sep="\t")
+            trdf.to_csv("treatments.csv", mode="a", index=False, sep="\t")
 
         if isinstance(cna_annotations, pd.DataFrame):
-            cna_annotations.to_csv("cna_annotated_cgi.csv", index=False, sep="\t", columns=['patient_id', 'sample_id', 'alteration', 'hugoSymbol', 'tumorType', 'oncogenic', 'mutationEffectDescription', 'gene_role', 'citationPMids', 'level_of_evidence', 'geneSummary', 'variantSummary', 'tumorTypeSummary'])
+            cna_annotations.to_csv(output, index=False, sep="\t", columns=['patient_id', 'sample_id', 'alteration', 'hugoSymbol', 'tumorType', 'oncogenic', 'mutationEffectDescription', 'gene_role', 'citationPMids', 'level_of_evidence', 'geneSummary', 'variantSummary', 'tumorTypeSummary'])
             trdf = pd.DataFrame(treatments)
-            trdf.to_csv("treatments_cgi_cna.csv", index=False, sep="\t")
+            trdf.to_csv("treatments.csv", mode="a", index=False, sep="\t")
 
         return 1
     else:
@@ -303,7 +301,7 @@ def generate_cgi_cna_file_from_list(genelist):
             file2.write(row)
         file2.close()
 
-def generate_temp_cgi_query_files(snv_annotations: pd.DataFrame = None, cna_annotations: pd.DataFrame = None, translocs: pd.DataFrame = None):
+def generate_temp_cgi_query_files(snv_annotations: pd.DataFrame = None, cna_annotations: pd.DataFrame = None, translocs: pd.DataFrame = None, append_to_annotations: bool = True):
     """
         Generate temporary CGI query files from annotations.
 
@@ -315,16 +313,27 @@ def generate_temp_cgi_query_files(snv_annotations: pd.DataFrame = None, cna_anno
     header = "chr\tpos\tref\talt\tsample\n"
     try:
         if isinstance(snv_annotations, pd.DataFrame):
-            with open("./tmp/snvs.ext", "w") as file1:
-                file1.write(header)
-
-                uniques = snv_annotations[['alteration']].drop_duplicates()
-                for indx, snv in uniques.iterrows():
-                    id = "SNV:"+snv['alteration']
-                    alt_split = snv['alteration'].split(':')
-                    row = alt_split[1]+'\t'+alt_split[2]+'\t'+alt_split[3]+'\t'+alt_split[4]+'\t'+id+'\n'
-                    file1.write(row)
-                file1.close()
+            if append_to_annotations:
+                with open("./tmp/snvs.ext", "w") as file1:
+                    file1.write(header)
+
+                    uniques = snv_annotations[['alteration']].drop_duplicates()
+                    for indx, snv in uniques.iterrows():
+                        id = "SNV:"+snv['alteration']
+                        alt_split = snv['alteration'].split(':')
+                        row = alt_split[1]+'\t'+alt_split[2]+'\t'+alt_split[3]+'\t'+alt_split[4]+'\t'+id+'\n'
+                        file1.write(row)
+                    file1.close()
+            else:
+                with open("./tmp/snvs.ext", "w") as file1:
+                    file1.write(header)
+
+                    uniques = snv_annotations[['hugoSymbol', 'chromosome', 'position', 'reference_allele', 'sample_allele', 'tumorType', 'referenceGenome']].drop_duplicates()
+                    for indx, snv in uniques.iterrows():
+                        id = "SNV:"+snv['hugoSymbol']+':'+snv['chromosome']+':'+str(snv['position'])+':'+snv['reference_allele']+':'+snv['sample_allele']
+                        row = snv['chromosome']+'\t'+str(snv['position'])+'\t'+snv['reference_allele']+'\t'+snv['sample_allele']+'\t'+id+'\n' #+'\t'+cryptocode.encrypt(snv.samples, settings.CRYPTOCODE)+'\n'
+                        file1.write(row)
+                    file1.close()
 
         if isinstance(cna_annotations, pd.DataFrame):
             header = "gene\tcna\tsample\n"

external_annotator.py

@@ -15,19 +15,21 @@
   --cgijobid <str>             Download results from CGI by jobid and apply annotations.
   --copy_number_alterations <str> Path to copy number alterations file.
   --somatic_variants <str>     Path to somatic variants file.
+  --output <str>               Path to output file.
+
 
 Examples:
-  python external_annotator.py --oncokbcna --copy_number_alterations path/to/cnas.tsv
-  python external_annotator.py --oncokbsnv --somatic_variants path/to/snvs.tsv
-  python external_annotator.py --cgiquery --somatic_variants path/to/snvs.tsv
-  python external_annotator.py --cgiquery --copy_number_alterations path/to/cnas.tsv
-  python external_annotator.py --cgiquery --cgijobid <jobid> --somatic_variants path/to/snvs.tsv
-  python external_annotator.py --cgiquery --cgijobid <jobid> --copy_number_alterations path/to/cnas.tsv
+  python external_annotator.py --oncokbcna --copy_number_alterations path/to/cnas.tsv --output path/to/output
+  python external_annotator.py --oncokbsnv --somatic_variants path/to/snvs.tsv --output path/to/output
+  python external_annotator.py --cgiquery --somatic_variants path/to/snvs.tsv --output path/to/output
+  python external_annotator.py --cgiquery --copy_number_alterations path/to/cnas.tsv --output path/to/output
+  python external_annotator.py --cgiquery --cgijobid <jobid> --somatic_variants path/to/snvs.tsv --output path/to/output
+  python external_annotator.py --cgiquery --cgijobid <jobid> --copy_number_alterations path/to/cnas.tsv --output path/to/output
 '''
 
 def main(**kwargs):
 
-
+    output = kwargs.get("output", ".")
     if kwargs["oncokbcna"] and kwargs["copy_number_alterations"]:
 
         cnas = pd.read_csv(kwargs["copy_number_alterations"], sep="\t")
@@ -46,7 +48,7 @@ def main(**kwargs):
         i = 0
         for c in chunks:
             i += 1
-            query_oncokb_cnas_to_csv(c, i)
+            query_oncokb_cnas_to_csv(c, output, i)
 
 
     if kwargs["oncokbsnv"] and kwargs["somatic_variants"]:
@@ -68,18 +70,18 @@ def main(**kwargs):
         i = 0
         for c in chunks:
             i += 1
-            query_oncokb_somatic_mutations(c, i)
+            query_oncokb_somatic_mutations(c, output, i)
 
     if kwargs["cgiquery"] and kwargs["somatic_variants"]:
         snvs = pd.read_csv(kwargs["somatic_variants"], sep="\t", dtype='string')
 
         if kwargs["cgijobid"]:
             jobid = kwargs["cgijobid"]
         else:
-            generate_temp_cgi_query_files(snvs, None, None)
-            jobid = launch_cgi_job_with_mulitple_variant_types("./tmp/snvs.ext",None, None, "OVSE", "hg38").replace('"', '')
+            generate_temp_cgi_query_files(snv_annotations=snvs)
+            jobid = launch_cgi_job_with_mulitple_variant_types(mutations_file="./tmp/snvs.ext", cancer_type="OVSE", reference="hg38").replace('"', '')
         time.sleep(30)
-        while query_cgi_job(jobid, snvs) == 0:
+        while query_cgi_job(jobid, output, snv_annotations=snvs) == 0:
             print("Waiting 30 seconds for the next try...")
             time.sleep(30)
 
@@ -89,11 +91,11 @@ def main(**kwargs):
         if kwargs["cgijobid"]:
             jobid = kwargs["cgijobid"]
         else:
-            generate_temp_cgi_query_files(None, cnas, None)
-            jobid = launch_cgi_job_with_mulitple_variant_types(None, "./tmp/cnas.ext", None, "OVSE", "hg38").replace('"', '')
+            generate_temp_cgi_query_files(cna_annotations=cnas)
+            jobid = launch_cgi_job_with_mulitple_variant_types(cnas_file="./tmp/cnas.ext", cancer_type="OVSE", reference="hg38").replace('"', '')
 
         time.sleep(30)
-        while query_cgi_job(jobid, None, cnas) == 0:
+        while query_cgi_job(jobid, output, cna_annotations=cnas) == 0:
             print("Waiting 30 seconds for the next try...")
             time.sleep(30)
 
@@ -106,6 +108,8 @@ def add_arguments(parser):
         parser.add_argument('--cgijobid', type=str, help='Download results from CGI by jobid')
         parser.add_argument('--copy_number_alterations', type=str, help='Path to copy number alterations file')
         parser.add_argument('--somatic_variants', type=str, help='Path to somatic variants file')
+        parser.add_argument('--output', type=str, default=".", help='Path to output directory for annotated files')
+
 
 
     parser = argparse.ArgumentParser()

local_annotator.py

@@ -135,9 +135,9 @@ def process_sample(gname, cna_grp, snv_grps, output, cores, ascats):
 
 if __name__ == "__main__":
     parser = argparse.ArgumentParser(description="Annotate genomic alterations.")
-    parser.add_argument("--output", type=str, required=True, help="Path to output file")
-    parser.add_argument("--somatic_variants", type=str, help="Path to somatic variants file")
-    parser.add_argument("--copy_number_alterations", type=str, help="Path to copy number alterations file")
+    parser.add_argument("--output", type=str, required=True, help="Path to output files")
+    parser.add_argument("--somatic_variants", type=str, help="Path to somatic variants files")
+    parser.add_argument("--copy_number_alterations", type=str, help="Path to copy number alterations files")
     parser.add_argument("--ascatestimates", type=str, required=True, help="Path to ASCAT estimates file")
     parser.add_argument("--cn_annotations", type=str, help="Path to filtered and annotated CNAs")
     parser.add_argument("--tumortype", type=str, default="HGSOC", help="Tumor type identifier (default: HGSOC)")

main.py

@@ -56,7 +56,6 @@ def run_external_annotator(args):
 
 
 def main():
-    parser = argparse.ArgumentParser(description="Run local and external annotators")
     parser = argparse.ArgumentParser(description="Run local and external annotators")
     parser.add_argument("--annotator", choices=["local", "external", "both"], required=True,
                         help="Choose which annotator to run")

oncokb_annotator.py

@@ -39,6 +39,7 @@ def handle_treatments_oncokb(jsondata, alt_type, alteration):
             'description': description,
             'treatment': drugs,
             'level_of_evidence': level,
+            'cgi_level':"",
             'citations': pmids,
             'tumorType': tumortype
         }))
@@ -65,7 +66,7 @@ def handle_drugs_field(jsondata):
         return None
 
 
-def query_oncokb_cnas_to_csv(cna_annotations: pd.DataFrame, i):
+def query_oncokb_cnas_to_csv(cna_annotations: pd.DataFrame, output, i):
 
     """
     Query OncoKB API to get annotations for copy number alterations (CNAs) and save the results to a CSV file.
@@ -146,16 +147,16 @@ def query_oncokb_cnas_to_csv(cna_annotations: pd.DataFrame, i):
             #print("Updated "+str(updatedf.count())+" CNAs")
         #cna_annotations.drop(columns=cna_annotations.columns[0], axis=1, inplace=True)
         header = False if i > 1 else True
-        cna_annotations.to_csv("cna_annotated_oncokb.csv", mode="a", index=False, header=header, sep="\t", columns=['patient_id', 'sample_id', 'alteration', 'hugoSymbol', 'tumorType', 'consequence', 'oncogenic', 'mutationEffectDescription', 'gene_role', 'citationPMids', 'level_of_evidence', 'geneSummary', 'variantSummary', 'tumorTypeSummary'])
+        cna_annotations.to_csv(output, mode="a", index=False, header=header, sep="\t", columns=['patient_id', 'sample_id', 'alteration', 'hugoSymbol', 'tumorType', 'consequence', 'oncogenic', 'mutationEffectDescription', 'gene_role', 'citationPMids', 'level_of_evidence', 'cgi_level', 'geneSummary', 'variantSummary', 'tumorTypeSummary'])
         trdf = pd.DataFrame(treatments)
-        trdf.to_csv("treatments_oncokb.csv", mode="a", header=header, index=False, sep="\t")
+        trdf.to_csv("treatments.csv", mode="a", header=header, index=False, sep="\t")
     else:
         print("Unable to request. Response: ", response.text)
 
     return response
 
 
-def query_oncokb_somatic_mutations(snv_annotations: pd.DataFrame, i):
+def query_oncokb_somatic_mutations(snv_annotations: pd.DataFrame, output, i):
     """
     Query OncoKB API to get annotations for somatic mutations and save the results to a CSV file.
 
@@ -225,9 +226,9 @@ def query_oncokb_somatic_mutations(snv_annotations: pd.DataFrame, i):
 
         print(snv_annotations)
         header = False if i > 1 else True
-        snv_annotations.to_csv("snv_annotated_oncokb.csv", mode="a", header=header, index=False, sep="\t", columns=['patient_id', 'sample_id', 'alteration', 'hugoSymbol', 'tumorType', 'consequence', 'oncogenic', 'mutationEffectDescription', 'gene_role', 'citationPMids', 'level_of_evidence', 'geneSummary', 'variantSummary', 'tumorTypeSummary'])
+        snv_annotations.to_csv(output, mode="a", header=header, index=False, sep="\t", columns=['patient_id', 'sample_id', 'alteration', 'hugoSymbol', 'tumorType', 'consequence', 'oncogenic', 'mutationEffectDescription', 'gene_role', 'citationPMids', 'level_of_evidence', 'cgi_level', 'geneSummary', 'variantSummary', 'tumorTypeSummary'])
         trdf = pd.DataFrame(treatments)
-        trdf.to_csv("treatments_oncokb_snv.csv", header=header, mode="a", index=False, sep="\t")
+        trdf.to_csv("treatments.csv", header=header, mode="a", index=False, sep="\t")
         #print("Updated " + str(len(snvdf)) + " CNAs")
     else:
         print("[ERROR] Unable to request. Response: ", print(response.text))

slurm_scripts/cna_annotation.sbatch

@@ -5,9 +5,6 @@
 #SBATCH --error=./logs/%A_%a.err
 #SBATCH --cpus-per-task=1
 #SBATCH --mem-per-cpu=5G
-#SBATCH	--partition=general,evmbig
-#SBATCH --exclude=evm06,evm07,evm08,evm09,evm10,evmfull01,evmbig
-#SBATCH --array=1-850
 
 echo $1
 ODANNOTATOR_PATH=""

slurm_scripts/snv_annotation.sbatch

@@ -5,9 +5,6 @@
 #SBATCH --error=./logs/%A_%a.err
 #SBATCH --cpus-per-task=1
 #SBATCH --mem-per-cpu=5G
-#SBATCH	--partition=general,evmbig
-#SBATCH --exclude=evm06,evm07,evm08,evm09,evm10,evmfull01,evmbig
-#SBATCH --array=1-850
 
 echo $1
 ODANNOTATOR_PATH=""