<jupyter_start><jupyter_text>Read Parsnp output<jupyter_code>with open("parsnp.xmfa", "r") as f:
    readlist = f.readlines()
def generate_read_list(reads):
    name_list = []
    read_list = []
    cluster = {}
    tempstr = ""
    for i in range(len(reads)):
        if reads[i][0] == "#":
            continue
        if reads[i][0] == ">":
            name_list.append(reads[i][:-1])
            # print(reads[i])
        elif i == len(reads) - 1:
            tempstr = tempstr+reads[i][:-1]
            read_list.append(tempstr)
            tempstr = ""
        elif reads[i+1][0] == ">":
            tempstr = tempstr+reads[i][:-1]
            read_list.append(tempstr)
            tempstr = ""
        else:
            tempstr = tempstr+reads[i][:-1]
    return (name_list, read_list)
(nlist, rlist) = generate_read_list(readlist)
clusters = {}
for i in range(len(nlist)):
    cluster = nlist[i].split(" ")[2]
    if cluster not in clusters:
        #exclude all LCBs
        if ('A' in rlist[i]) or ('T' in rlist[i]) or ('C' in rlist[i]) or ('G' in rlist[i]) or ('-' in rlist[i]):
            continue
        clusters.update({cluster: rlist[i]})
os.system("mkdir ./All_MUMs")
for i in clusters:
    with open("./All_MUMs/" + i + ".fasta", "a+") as f:
        f.write(">"+i+"\n")
        f.write(clusters[i])
clusterl = os.listdir("./All_MUMs/")<jupyter_output><empty_output><jupyter_text>Run Blast for all MUMs ** Important: Some blastn result can be None, delete those files manually.<jupyter_code>for i in tq(clusterl):
    os.system("blastn -max_target_seqs 1000 -db /home/Users/yf20/ncbi_database/newnt/nt -query "+ "./All_MUMs/" + i + " -out " + i + ".out -outfmt '6 qseqid sseqid pident evalue stitle' -num_threads 70")
os.system("mkdir ./newblastresult/")
os.system("mv *.out ./newblastresult/")
blastresultfile = os.listdir("newblastresult/")
strain = []
percent2strain = []
for i in blastresultfile:
    with open ("newblastresult/" + i, "r") as f:
        l = f.readlines()
        content = []
        p2s = {}
        for j in l:
            k = j.split("\t")
            content.append(k)
            stranname = k[-1][:-1]
            if stranname not in strain:
                strain.append(stranname)
            strpercent = k[2]
            if stranname not in p2s:
                p2s.update({stranname: float(strpercent)})
            read = clusters[k[0]]
        p2s.update({"MUM":read})
        percent2strain.append((p2s))<jupyter_output><empty_output><jupyter_text>Group different, important(required) Strains<jupyter_code>mtbs = []
unknowns = []
mbovis = []
mnontbs = []
afs = []
for i in strain:
    si = i.split(" ")
    if si[0] == "Mycobacterium" and si[1] == "tuberculosis":
        #Mtb Bovis strain is also bovis strain
        if si[2] == "variant" and si[3] == "bovis":
            mbovis.append(i)
        else:
            mtbs.append(i)
    elif si[1] == "sp.":
        unknowns.append(i)
    elif si[0] == "Mycobacterium" and si[1] == "bovis":
        mbovis.append(i)
    elif si[0] == "Mycobacterium" and si[1] == "africanum":
        afs.append(i)
    else:
        mnontbs.append(i)
whole = mtbs+mbovis+afs+mnontbs+unknowns
res = pd.DataFrame(columns=["MUM"]+whole)
for i in percent2strain:
    res = res.append(i, ignore_index=True)
res = res.fillna(0)
res.set_index(["MUM"], inplace = True)<jupyter_output><empty_output><jupyter_text>Get all 16S strains<jupyter_code>Strain16S = res
for i in strain:
    keys = i.split(" ")
    if "16S" in keys:
        continue
    Strain16S = Strain16S.drop(i, axis=1)<jupyter_output><empty_output><jupyter_text>See All Complete Strains<jupyter_code>rescg = res
for i in strain:
    keys = i.split(" ")
    if keys[-1] == "genome" and keys[-2] == "complete" and keys[0] == "Mycobacterium" and keys[1] != "sp.":
        continue
    rescg = rescg.drop(i, axis=1)<jupyter_output><empty_output><jupyter_text>See Bovis<jupyter_code>boviscg = []
for i in mbovis:
    keys = i.split(" ")
    if keys[-1] == "genome" and keys[-2] == "complete" and keys[0] == "Mycobacterium" and keys[1] != "sp.":
        boviscg.append(i)
resbovis = rescg[boviscg]
tbcg = []
for i in mtbs:
    keys = i.split(" ")
    if keys[-1] == "genome" and keys[-2] == "complete" and keys[0] == "Mycobacterium" and keys[1] != "sp.":
        tbcg.append(i)
resnontb = rescg.drop(tbcg, axis=1)<jupyter_output><empty_output><jupyter_text>See All Unique to TB<jupyter_code>MUMs = list(resnontb.index.values)
droplist = []
for i in MUMs:
    for j in (resnontb.loc[i]):
        if j == 100 and (i not in droplist):
            droplist.append(i)
uniquttotb = rescg.drop(droplist)
with open("all_strains.csv", "w") as f:
    f.write(res.to_csv())
with open("complete_genomes.csv", "w") as f:
    f.write(rescg.to_csv())
with open("bovis_complete_genomes_strains.csv", "w") as f:
    f.write(resbovis.to_csv())
with open("non_tb_complete_genomes_strains.csv", "w") as f:
    f.write(resnontb.to_csv())
with open("unique_to_tb.csv", "w") as f:
    f.write(uniquttotb.to_csv())
with open("16S.csv", "w") as f:
    f.write(Strain16S.to_csv())<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.1666015625, 'fraction_numerical': 0.0142578125, 'mean_word_length': 2.9271472392638036, 'pattern_counts': {'":': 12, '<': 29, '<?xml version=': 0, '>': 32, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '7', 'dir': 'jupyter-structured-clean-dedup', 'id': '412142', 'license': 'MIT', 'n_tokens_mistral': 2024, 'n_tokens_neox': 1933, 'n_words': 415, 'path': 'Notebook_Templates/resultparse.ipynb', 'repo_name': 'Fu-Yilei/SUEDE', 'size': '5120'}