Commit 07b6f626 authored by Roman Sarrazin-Gendron's avatar Roman Sarrazin-Gendron
Browse files

fixed rotations and inline ss

parent a68b7826
......@@ -480,13 +480,14 @@ def parse_sequence(seq, modules, ss, dataset, BNs, t=-10, samplesize=20000, Lamb
#print(seq)
#print(struct)
if not str(r.node.positions) in modules_predicted[r.module.ID]:
#print("about to score",r.seq, r.rotation, r.pos)
r.eval(fold_compound=fc)
#print("SCORED MODULE",r.module.ID, r.seq, r.pos, r.score)
#print("PASSING THRESHOLD",r.score,t)
if r.score > t:
#print(r.module.ID, r.seq, r.pos, r.score)
#print("Adding module match to results")
modules_predicted[r.module.ID][str(r.node.positions)] = [r.seq, r.pos, prob,r.node.positions, r.score]
modules_predicted[r.module.ID][str(r.node.positions)] = [r.seq, r.pos, prob,r.node.positions, r.score, r.mapping]
else:
#print("ADDING PREDICTION TO OTHERS",r.pos, modules_predicted[r.module.ID])
......
......@@ -743,7 +743,7 @@ if __name__ == "__main__":
#print("THE HITS")
#print(all_results[inputSeqKey])
modules_in_svg, chef_ss = bp_chefs_choice(all_results[inputSeqKey],seqInfo[seqCounter],arguments["t"],finalName)
modules_in_svg, chef_ss = bp_chefs_choice(all_results[inputSeqKey],seqInfo[seqCounter],float(arguments["t"]),finalName)
#now we need to fill svg hits
......
......@@ -325,7 +325,9 @@ class RNAModule(BayesianModel):
def eval(self, sequence):
#print('testing seq against nodes',sequence,self.nodes())
candidate = self.candidate_from_sequence(sequence)
#print('testing seq against nodes',sequence,self.nodes(),candidate)
if candidate is None:
return [-100, None]
# in the bayesnet object, the nucleotide are stored as number
......@@ -342,6 +344,7 @@ class RNAModule(BayesianModel):
this_node_evidence_dict = {}
for anc in ancestors_of_node:
this_node_evidence_dict[anc] = candidate[anc]
#print("evaluating",node, candidate[node], this_node_evidence_dict)
node_probability = math.log(compute_probability(self, node, candidate[node], this_node_evidence_dict))
node_probabilities.append(node_probability)
module_probability = module_probability + node_probability
......@@ -400,6 +403,22 @@ class RNAModule(BayesianModel):
return candidate
def mapping_from_sequence(self, sequence):
#print("SCORING CANDIDATE", self.ID, self.nodes(),sequence, "GAPS",self.gaps)
#takes as input a sequence with stars, and outputs a dictionary of node values
candidate = {}
if "T" in sequence:
sequence.replace("T","U")
sequence = list(sequence)
sequence = [i for i in sequence if i.isalpha()]
for ind,node in enumerate(self.nodes(data=False)):
candidate[node] = sequence[ind]
return candidate
def eval_constraint_folding(self,seq_score,full_seq,positions,rotated_ss, rotated_gaps, fc):
#print('testing candidate with constraint folding',rotated_ss, rotated_ss.split("*"))
......
......@@ -58,7 +58,7 @@ class ExactScanResult:
if rotation:
self.rotate(rotation)
#print("LEN POS",self.pos,"LEN SEQ",self.seq, "ROTATION",rotation)
#print("DONE ROTATION, CURRENT SEQUENCE", self.seq)
#TODO : need to think about adjusting POSITIONS to gaps
self.gaps = self.module.gaps_per_strand
......@@ -66,7 +66,8 @@ class ExactScanResult:
#print("PREFIX SEQUENCES",self.seq)
self.reintegrate_module_gaps_to_sse_sequence_and_positions()
#print("Reintegrated module graps, current seq", self.seq)
self.mapping = {}
def rotate_stackings(self):
......@@ -112,7 +113,8 @@ class ExactScanResult:
if ind not in gaps[0]:
new_seq += nuc
else:
for strand,subseq in enumerate(self.rotate_list(seq.split("&"))):
#for strand,subseq in enumerate(self.rotate_list(seq.split("&"))):
for strand,subseq in enumerate(seq.split("&")):
#print("STRAND",strand,"GAPS",gaps, seq, seq.split("&"))
if not gaps[strand]:
new_seq+=subseq
......@@ -147,7 +149,7 @@ class ExactScanResult:
"""
#TODO: need to make sure the gaps are managed when scoring the sequence, specifically through rotation
#print("info",self.module.stackings,self.module.ID,self.module.n_nodes,self.seq)
#print("info eval",self.rotation,self.pos,self.module.ID,self.module.n_nodes,self.seq)
try:
return self.score
except:
......@@ -161,6 +163,7 @@ class ExactScanResult:
return -1
# print(self.module.ID,self.node,self.pos)
score = self.module.eval(self.seq)[0]
self.mapping = self.module.mapping_from_sequence(self.seq)
#print("COMPUTING SEQUENCE SCORE",score, self.pos, self.seq)
self.score = score
......@@ -182,6 +185,8 @@ class ExactScanResult:
else:
#print("Candidate scoring", tmp_seq,self.module.eval(tmp_seq)[0] )
total_score+= self.module.eval(tmp_seq)[0]
if len(self.mapping)==0:
self.mapping = self.module.mapping_from_sequence(tmp_seq)
mean_score = total_score/(len(sequences)-gapped_cand)
# print("score for this alignment:",mean_score)
self.score = mean_score
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment