猪的stop数据库芯片构建

import pandas as pd

name_lst = []
allseq_lst = []
seq_lst = ["test"]
with open("AllExonSeq.fa", "r") as f:
    for line in f:
        if line.startswith(">"):
            allseq_lst.append("".join(seq_lst))
            seq_lst = []
            name = line.strip("\n")
            name_lst.append(name)
        else:
            seq_lst.append(line.strip("\n"))
allseq_lst.remove("test")
allseq_lst.append(seq_lst)
df = pd.DataFrame({"seqName": name_lst, "sequence": allseq_lst})
df2 = pd.read_csv("AllExonInfoTop3.txt", sep="\t", header=0)
df = df.set_index("seqName")
df2 = df2.set_index("seqName")
df3 = df2.join(df, how="inner")
df3 = df3.reset_index()
df3.drop_duplicates("seqName", keep="first", inplace=True)
df3.to_csv("exonWithAnno.txt", sep="\t", index=None)
ss = df3["seqName"].tolist()
sn = df3["sequence"].tolist()
with open("AllGeneTop3Exon.fa", "a") as l:
    for name, seq in zip(ss, sn):
        l.write("{}\n{}\n".format(name, seq))

import pandas as pd


def complement(seq):
    return seq.translate(str.maketrans("ACGT", "TGCA"))


def revcomp(seq):
    return complement(seq)[::-1]


df = pd.read_csv("exonWithAnno.txt", sep="\t", header=0)
df2 = df[df["Strand"] < 0]
df1 = df[df["Strand"] > 0]
df2["sequence"] = df2["sequence"].apply(revcomp)
fin = pd.concat([df1, df2], axis=0, ignore_index=True)
ss = fin["seqName"].tolist()
nn = fin["sequence"].tolist()
with open("allexoninzhengliantop3.fa", "a") as f:
    for name, seq in zip(ss, nn):
        f.write(f"{name}\n{seq}\n")

import pandas as pd
import sys
name_header = sys.argv[1]
gRNA_filer = pd.read_csv("./tmp/{}.first_filter.txt".format(name_header), sep="\t", header=0)
gRNA_Z = gRNA_filer[((gRNA_filer["orientation"] == "FWD") & (gRNA_filer["strand"] == "+")) | ((gRNA_filer["orientation"] == "RVS") & (gRNA_filer["strand"] == "-"))].reset_index()[['gene','target']]
gRNA_F = gRNA_filer[((gRNA_filer["orientation"] == "RVS")&(gRNA_filer["strand"] == "+")) | ((gRNA_filer["orientation"] == "FWD")&(gRNA_filer["strand"] == "-"))].reset_index()[['gene','target']]
gRNA_Z.to_csv("./tmp/{}.part_1.txt".format(name_header), sep="\t", index=None, header=None)
gRNA_F.to_csv("./tmp/{}.part_2.txt".format(name_header), sep="\t", index=None, header=None)

import re
import sys
name_header = sys.argv[1]
def read_cds(fs):
    name_lst = []
    seq_lst=[]
    with open(fs, 'r') as f:
        for line in f:
            if line.startswith(">"):
                name = line.strip("\n").replace(">", "")
                name_lst.append(name)
            else:
                seq_lst.append(line.strip("\n"))
    fin=[(key,value) for key,value in zip(name_lst,seq_lst) if 'unavailable' not in value]
    return fin

def complement(seq):
    return seq.translate(str.maketrans('ACGT', 'TGCA'))

def revcomp(seq):
    return complement(seq)[::-1]

def motif_Z(gRNA, site):
    if gRNA[:20][site: site + 3] == "CAG":
        return "CAG"
    elif gRNA[:20][site: site + 3] == "CGA":
        return "CGA"
    elif gRNA[:20][site: site + 3] == "CAA":
        return "CAA"
    else:
        return 0
def stop_test_for_Z(n, cds_lst, gRNA):
    new_cds_lst = [cds for cds in cds_lst if gRNA in cds]
    score_lst = []
    if len(new_cds_lst) > 0:
        for haha in new_cds_lst:
            m = haha.find(gRNA) + 1
            score = (m - 1 + n - 1)
            score_lst.append(score)
        if len([new_score for new_score in score_lst if score % 3 == 0]) > 0:
            return 1
        else:
            return 0
    else:
        return 0
def motif_F(gRNA, site):
    if gRNA[:20][site: site + 3] == "CCA":
        return 1
    elif gRNA[:20][site-1: site + 2] == "CCA":
        return 2
    else:
        return 0
def stop_test_for_F(n, cds_lst, gRNA):
    gRNA = revcomp(gRNA)
    new_cds_lst = [cds for cds in cds_lst if gRNA in cds]
    score_lst = []
    if len(new_cds_lst) > 0:
        for haha in new_cds_lst:
            m = haha.find(gRNA) + 1 - 1 + 23
            score = (m - n)
            score_lst.append(score)
        if len([new_score for new_score in score_lst if score % 3 == 0]) > 0:
            return 1
        else:
            return 0
    else:
        return 0
def detective_Z(gRNA, Gene):
    #cds_lst=all_cds[all_cds['Gene']==Gene]['seq'].tolist()
    cds_lst=[value for (key,value) in all_cds if key==Gene]
    stopSitelst=[]
    for m in re.finditer("C", gRNA[:20]):
        if m.start() + 1 >= 4 and m.start() + 1 <= 8:
            if motif_Z(gRNA, m.start()):
                if stop_test_for_Z((m.start() + 1), cds_lst, gRNA):
                    stopSitelst.append(str(m.start() + 1))
    if len(stopSitelst) > 0:
        fin="/".join(stopSitelst)
    else:
        fin="None"
    return fin

def detective_F(gRNA, Gene):
    #cds_lst=all_cds[all_cds['Gene']==Gene]['seq'].tolist()
    cds_lst=[value for (key,value) in all_cds if key==Gene]
    stopSitelst=[]
    for m in re.finditer("C", gRNA[:20]):
        if m.start() + 1 >= 4 and m.start() + 1 <= 8:
            if motif_F(gRNA, m.start()):
                if stop_test_for_F(((gRNA[:20][m.start() - 1 : m.start() + 3]).find("CCA") + 1 + 2), cds_lst, gRNA):
                    stopSitelst.append(str(m.start() + 1))
    if len(stopSitelst) > 0:
        fin="/".join(stopSitelst)
    else:
        fin="None"
    return fin

all_cds = read_cds("AllGeneCDS2.fa")
shunshi=[]
with open("./tmp/{}.part_1.txt".format(name_header),'r') as f:
    for line in f:
        gene,target=line.split("\t")
        target=target.strip("\n")
        shunshi.append((gene,target))

new_lst=[(key,value,detective_Z(value,key)) for (key,value) in shunshi]

with open("./tmp/{}.shunshi.out".format(name_header),"a") as f:
    for (gene,target,site) in new_lst:
        f.write("{}\t{}\t{}\n".format(gene,target,site))

fanshi=[]
with open("./tmp/{}.part_2.txt".format(name_header),'r') as f:
    for line in f:
        gene,target=line.split("\t")
        target=target.strip("\n")
        fanshi.append((gene,target))

new_lst2=[(key,value,detective_F(value,key)) for (key,value) in fanshi]

with open("./tmp/{}.fanshi.out".format(name_header),"a") as f:
    for (gene,target,site) in new_lst2:
        f.write("{}\t{}\t{}\n".format(gene,target,site))

import pandas as pd
import sys
name_header = sys.argv[1]
df1 = pd.read_csv("./tmp/{}.shunshi.out".format(name_header), sep="\t", header=None)
df2 = pd.read_csv("./tmp/{}.fanshi.out".format(name_header), sep="\t", header=None)
fin = pd.concat([df2, df1], ignore_index=True)
fin.columns = ['gene', 'target', 'site']
fin2 = fin[fin['site'] != "None"]
fin2.to_csv("./tmp/{}.second_filter.txt".format(name_header),sep="\t",index=None)