Suggested Bestimate fixes,

.gitignore

@@ -2,3 +2,4 @@
 .Rhistory
 .RData
 .Ruserdata
+scratch

BEstimate/BEstimate.py

@@ -12,7 +12,7 @@
 from Bio import SeqIO
 from Bio import pairwise2
 from Bio.pairwise2 import format_alignment
-
+import numpy as np
 
 # -----------------------------------------------------------------------------------------#
 # Take inputs
@@ -2386,7 +2386,7 @@ def summarise_guides(last_df):
 		summary_df.loc[i, "gRNA_Target_Sequence"] = ";".join(
 			[x for x in list(guide_df.gRNA_Target_Sequence.unique()) if x is not None])
 		summary_df.loc[i, "gRNA_flanking_sequences"] = ";".join(
-			[x for x in list(guide_df.gRNA_flanking_sequences.unique()) if x is not None])
+			[str(x) for x in list(guide_df.gRNA_flanking_sequences.unique()) if x is not None and isinstance(x, float) and not np.isnan(x)])
 		summary_df.loc[i, "gRNA_Target_Location"] = ";".join(
 			[x for x in list(guide_df.gRNA_Target_Location.unique()) if x is not None])
 		summary_df.loc[i, "Edit_Location"] = ";".join(

BEstimate/add_flaking.py

@@ -9,7 +9,8 @@
 
 # Import necessary packages
 import os, sys, pandas, argparse, requests
-
+import time
+import itertools
 
 # -----------------------------------------------------------------------------------------#
 # Take inputs
@@ -54,7 +55,7 @@ def take_input():
 # -----------------------------------------------------------------------------------------#
 # Functions
 
-def extract_gRNA_flan_sequence(location, direction, fivep, threep):
+def extract_gRNA_flan_sequence(location, direction, fivep, threep, assembly):
 	"""
 	Annotating flanking regions of the gRNAs
 	:param location: Location of the gRNA target sequence
@@ -68,29 +69,49 @@ def extract_gRNA_flan_sequence(location, direction, fivep, threep):
 		grna_strand = "-1"
 	else:
 		grna_strand = "1"
+
+	server = "http://grch37.rest.ensembl.org" if assembly == "hg19" else "https://rest.ensembl.org"
 
-	grna_flan_ensembl = self.server + "/sequence/region/human/%s:%s?expand_3prime=%s;expand_5prime=%s;content-type=text/plain" \
-						% (location, grna_strand, threep, fivep)
-	grna_flan_request = requests.get(grna_flan_ensembl,
-									 headers={"Content-Type": "text/plain"})
+	# X:154444204-154444226;X:154444207-154444229
+	locations = location.split(";")
+	resp = []
+	for _loc in locations:
 
-	if grna_flan_request.status_code != 200:
-		print("No response from ensembl sequence!\n")
-		return "API problem"
-	else:
-		return grna_flan_request.text
+		grna_flan_ensembl = server + "/sequence/region/human/%s:%s?expand_3prime=%s;expand_5prime=%s;content-type=text/plain" \
+							% (_loc, grna_strand, threep, fivep)
+		grna_flan_request = requests.get(grna_flan_ensembl,
+											headers={"Content-Type": "text/plain"})
+
+		while grna_flan_request.status_code == 429:
+			print("Received status 429. Waiting for 1 second before sending next request.")
+			time.sleep(1)
+			grna_flan_request = requests.get(grna_flan_ensembl, headers={"Content-Type": "text/plain"})
+
+		if grna_flan_request.status_code != 200:
+			print(f"No response from ensembl sequence! {grna_flan_request.status_code}\n")
+			resp.append("API problem")
+		else:
+			resp.append(grna_flan_request.text)
 
+	return ";".join(resp)
 
 
 def main(path, args):
 
 	df = pandas.read_csv(path + args["FILE"] + ".csv", index_col=0)
-
+	flan_5 = args["FLAN_5"]
+	flan_3 = args["FLAN_3"]
+	assembly = args["ASSEMBLY"]
 	df["gRNA_flanking_sequences"] = df.apply(
 		lambda x: extract_gRNA_flan_sequence(location=x.CRISPR_PAM_Location, direction=x.Direction,
-											 fivep=flan_5, threep=flan_3), axis=1)
-
-	if "API problem" in df.gRNA_flanking_sequences.unique():
+											 fivep=flan_5, threep=flan_3, assembly = assembly), axis=1)
+
+	unique_responses = list(
+		itertools.chain(
+			*map(lambda x: x.split(";"), df.gRNA_flanking_sequences.unique())
+		)
+	)
+	if "API problem" in unique_responses:
 		print("API related problems")
 		return False
 	else:
@@ -102,10 +123,9 @@ def main(path, args):
 if __name__ == '__main__':
 
 	args = take_input()
-
 	# Path
 	path = ""
 	if args["PATH"][-1] == "/": path = args["PATH"]
 	else: path = args["PATH"] + "/"
 
-	_ = main(args)
+	_ = main(path, args)

bestimate.yml

@@ -1,9 +1,7 @@
 name: bestimate
 channels:
   - bioconda
-  - anaconda
   - conda-forge
-  - defaults
 dependencies:
   - _libgcc_mutex=0.1=main
   - _openmp_mutex=5.1=1_gnu