<div dir="ltr"><div class="markdown-here-wrapper" style=""><p style="margin:0px 0px 1.2em!important">Hi Sean,</p>
<p style="margin:0px 0px 1.2em!important">If you are doing this using a for loop, won’t your <code style="font-size:0.85em;font-family:Consolas,Inconsolata,Courier,monospace;margin:0px 0.15em;padding:0px 0.3em;white-space:pre-wrap;border:1px solid rgb(234,234,234);border-radius:3px;display:inline;background-color:rgb(248,248,248)">SeqIO.write(new_rec,"NewSeq.fasta","fasta")</code> simply write the last record of the loop.<br>I believe you should be passing a handle to SeqIO.write, with an append flag.</p>
<p style="margin:0px 0px 1.2em!important">Also, since you are creating a new <code style="font-size:0.85em;font-family:Consolas,Inconsolata,Courier,monospace;margin:0px 0.15em;padding:0px 0.3em;white-space:pre-wrap;border:1px solid rgb(234,234,234);border-radius:3px;display:inline;background-color:rgb(248,248,248)">Seq</code> object from a new string altogether, I don’t think <code style="font-size:0.85em;font-family:Consolas,Inconsolata,Courier,monospace;margin:0px 0.15em;padding:0px 0.3em;white-space:pre-wrap;border:1px solid rgb(234,234,234);border-radius:3px;display:inline;background-color:rgb(248,248,248)">tomutable</code> is required, unless you want to do something like this(the `append_handle` initialization being out of scope of for loop):</p>
<pre style="font-size:0.85em;font-family:Consolas,Inconsolata,Courier,monospace;font-size:1em;line-height:1.2em;margin:1.2em 0px"><code style="font-size:0.85em;font-family:Consolas,Inconsolata,Courier,monospace;margin:0px 0.15em;padding:0px 0.3em;white-space:pre-wrap;border:1px solid rgb(234,234,234);border-radius:3px;display:inline;background-color:rgb(248,248,248);white-space:pre;overflow:auto;border-radius:3px;border:1px solid rgb(204,204,204);padding:0.5em 0.7em;display:block!important">append_handle = open('new_records.fasta','a')
mutable_rec = refseq.tomutable()
for index,value in enumerate(mutable_rec):
    mutable_rec.seq[i] = filter_coverage(index,value)
SeqIO.write(mutable_rec, append_handle, 'fasta')
</code></pre><div title="MDH:SGkgU2Vhbiw8YnI+PGJyPklmIHlvdSBhcmUgZG9pbmcgdGhpcyB1c2luZyBhIGZvciBsb29wLCB3
b24ndCB5b3VyIGBTZXFJTy53cml0ZShuZXdfcmVjLCJOZXdTZXEuPHdicj5mYXN0YSIsImZhc3Rh
IilgIHNpbXBseSB3cml0ZSB0aGUgbGFzdCByZWNvcmQgb2YgdGhlIGxvb3AuPGRpdj5JIGJlbGll
dmUgeW91IHNob3VsZCBiZSBwYXNzaW5nIGEgaGFuZGxlIHRvIFNlcUlPLndyaXRlLCB3aXRoIGFu
IGFwcGVuZCBmbGFnLjwvZGl2PjxkaXY+PGJyPjwvZGl2PjxkaXY+QWxzbywgc2luY2UgeW91IGFy
ZSBjcmVhdGluZyBhIG5ldyBgU2VxYCBvYmplY3QgZnJvbSBhIG5ldyBzdHJpbmcgYWx0b2dldGhl
ciwgSSBkb24ndCB0aGluayBgdG9tdXRhYmxlYCBpcyByZXF1aXJlZCwgdW5sZXNzIHlvdSB3YW50
IHRvIGRvIHNvbWV0aGluZyBsaWtlIHRoaXM6PC9kaXY+PGRpdj48YnI+PC9kaXY+PGRpdj5gYGA8
L2Rpdj48ZGl2PmFwcGVuZF9oYW5kbGUgPSBvcGVuKCduZXdfcmVjb3Jkcy5mYXN0YScsJ2EnKTwv
ZGl2PjxkaXY+bXV0YWJsZV9yZWMgPSByZWZzZXEudG9tdXRhYmxlKCk8L2Rpdj48ZGl2PmZvciBp
bmRleCx2YWx1ZSBpbiBlbnVtZXJhdGUobXV0YWJsZV9yZWMpOjwvZGl2PjxkaXY+Jm5ic3A7ICZu
YnNwOyBtdXRhYmxlX3JlY1tpXSA9IGZpbHRlcl9jb3ZlcmFnZShpbmRleCx2YWx1ZSk8L2Rpdj48
ZGl2PlNlcUlPLndyaXRlKG11dGFibGVfcmVjLCBhcHBlbmRfaGFuZGxlLCAnZmFzdGEnKTxicj48
ZGl2PjxkaXY+YGBgPC9kaXY+PC9kaXY+PC9kaXY+PGRpdj48YnI+PC9kaXY+" style="height:0;width:0;max-height:0;max-width:0;overflow:hidden;font-size:0em;padding:0;margin:0">​</div></div></div><div class="gmail_extra"><br><div class="gmail_quote">On 16 March 2016 at 13:36, Sean Brimer <span dir="ltr"><<a href="mailto:skbrimer@gmail.com" target="_blank">skbrimer@gmail.com</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div class="HOEnZb"><div class="h5"><div dir="ltr"><span style="font-size:12.8px">Hello Biopython users, </span><div style="font-size:12.8px"><br></div><div style="font-size:12.8px">I'm trying to create a script that will allow me to add 'N' to a final fasta of a template-guided assembly. I have already posted this question on BioStars and received some good feedback and with it I feel I'm 99% done with what I am trying to complete. </div><div style="font-size:12.8px"><br></div><div style="font-size:12.8px">I am using python3 and Biopython 1.66</div><div style="font-size:12.8px"><br></div><div style="font-size:12.8px">The link to the BioStars post is here: <a href="https://www.biostars.org/p/179755/" target="_blank">https://www.biostars.org/p/179755/</a></div><div style="font-size:12.8px"><br></div><div style="font-size:12.8px">The code seems to do fine on making a dictionary of the coverage file and then making a set containing only the "no coverage" areas. I can edit the fasta header information i.e. new id, new description however I am unable to make an editted version of the Seq object. I have tried using the str() function and the .tomutable() function in several ways but none of them seem to be giving me my expected output. I suspect the error is in how I'm trying to call/edit the Seq function but from the tutorial guide I'm not sure why it's not working. </div><div style="font-size:12.8px"><br></div><div style="font-size:12.8px">Any guidance would be appreciated.</div><div style="font-size:12.8px"><br></div><div style="font-size:12.8px">Thank you, </div><div style="font-size:12.8px">Sean</div><div style="font-size:12.8px"><br></div><div style="font-size:12.8px">my code is below:</div><blockquote class="gmail_quote" style="font-size:12.8px;margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex"><br>from Bio.SeqRecord import SeqRecord<br>from Bio.Alphabet import IUPAC<br>from Bio.Seq import Seq<br>from Bio import SeqIO<br>import argparse <br></blockquote><div style="font-size:12.8px"> </div><blockquote class="gmail_quote" style="font-size:12.8px;margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex">## this makes a dict of the samtools depth coverage input file of all the 0's<br>coverageDict = {}<br></blockquote><div style="font-size:12.8px"> </div><blockquote class="gmail_quote" style="font-size:12.8px;margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex">## this loops over the input depth information and appends the dictionary with <br>## a key,value as genome position, depth of coverage<br>coverage = open("/home/sbrimer/Desktop/Columbia ST/1819_1/coverage.txt","r")<br>for line in coverage:<br>    coverages = line.split("\t")<br>    coverageDict[coverages[1]] = coverages[2]<br>coverage.close()<br></blockquote><div style="font-size:12.8px"> </div><blockquote class="gmail_quote" style="font-size:12.8px;margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex">## this should read in the sequence file as an oblect, change it to a mutable<br>## object, then anywhere the coverageDict value = 0. Replace that nt with "N"<br></blockquote><div style="font-size:12.8px"> </div><blockquote class="gmail_quote" style="font-size:12.8px;margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex">missing = {index for index,value in coverageDict.items() if value == 0}<br>def filter_coverage(index,value):<br>    return 'N' if index in missing else value<br></blockquote><div style="font-size:12.8px"> </div><blockquote class="gmail_quote" style="font-size:12.8px;margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex">newSeqRec = []<br>with open("newfile.fasta","w") as handle:<br>    for seq_record in SeqIO.parse("Salmonella enterica subsp.fasta","fasta"):<br>        refseq = seq_record.seq<br>        newSeq = "".join(filter_coverage(index,value) for index,value in enumerate(refseq.tomutable(), start=1))<br>        newSeqRec = SeqRecord(Seq(newSeq),id ="new_id", description = str(seq_record.description))<br>    SeqIO.write(newSeqRec, handle, "fasta")<br>handle.close()</blockquote></div>
</div></div><br>_______________________________________________<br>
Biopython mailing list  -  <a href="mailto:Biopython@mailman.open-bio.org">Biopython@mailman.open-bio.org</a><br>
<a href="http://mailman.open-bio.org/mailman/listinfo/biopython" rel="noreferrer" target="_blank">http://mailman.open-bio.org/mailman/listinfo/biopython</a><br></blockquote></div><br></div>