在 NGS 数据中寻找 CDR答案

【问题标题】：finding CDRs in NGS data在 NGS 数据中寻找 CDR
【发布时间】：2023-02-06 14:40:08
【问题描述】：

我有数百万个 fasta 格式的序列，想提取 CDR（CDR1、CDR2 和 CDR3）。我只选择了一个序列作为示例，并尝试提取 CDR1，但无法提取 CDR1。

顺序：-'FYSHSAVTLDESGGGLQTPGGGLSLVCKASGFTFSSYGMMWVRQAPGKGLEYVAGIRNDA GDKRYGSAVQGRATISRDNGQSTVRLQLNNLRAEDTGTYFCAKESGCYWDSTHCIDAWGH GTEVIVSTGG'.

cdr1 开始于：-'VCKASGFTFS'，最多替换三个，但必须是第二名 C。cdr1 结束于：-'WVRQAP'，最多有两个替换，但必须在第 3 名的 R。

提取的cdr1应该是SYGMM

def cdr1_in(cdr_in): #VCKASGFTFS
    pin=0
    max_pin=3       
    
    if cdr[1]!='C':
        pin+=1
    if cdr[0]!='V':
        pin+=1
    if cdr[2]!='K':
        pin+=1
    if cdr[3]!='A':
        pin+=1    
    if cdr[4]!='S':
        pin+=1
    if cdr[5]!='G':
        pin+=1
    if cdr[6]!='F':
        pin+=1
    if cdr[7]!='T':
        pin+=1    
    if cdr[8]!='F':
        pin+=1
    if cdr[9]!='S':
        pin+=1   
  
    if pin<max_pin:
        print('CDR_in pattern', cdr_in)
        # print('CDR_starts from', arr.index(cdr_in)+9)
        return (arr.index(cdr_in)+9)
 
    def cdr1_out(cdr_out):#WVRQAP
    
        pin=0
        max_pin=2            
        if cdr[1]!='V':
            pin+=1
        if cdr[0]!='W':
            pin+=1
        if cdr[2]!='R':
            pin+=1
        if cdr[3]!='Q':
            pin+=1    
        if cdr[4]!='A':
            pin+=1
        if cdr[5]!='P':
            pin+=1
            
        if pin<max_pin:
            # print('CDR_in pattern', cdr_out)
            # print('CDR_ends at', arr.index(cdr_out))
            return (arr.index(cdr_out))
 

K=10
arr=sequence
for i in range(len(arr)-k+1):
        slider=arr[i:k+i]
        print("CDR_1 is:", arr[cdr1_in(slider): cdr1_out(slider)])

【问题讨论】：

标签： python regex python-re biopython fuzzy-search

【解决方案1】：

我假设您正在分析免疫测序数据并且 CDR 是指 B 或 T 细胞受体的互补决定区域是否正确？数据来自人类还是小鼠？如果是这种情况，与其重新发明轮子，不如看看现有的工具。我用过mixcr。另一个流行的工具是IMGT/HighV-QUEST，但据我所知，它只能作为网络应用程序使用，不能用于大型数据集。如果它们不符合您的目的，您至少可以得到有关如何进行的提示。

【讨论】：

谢谢你的建议。我已经使用 IMGT/KABAT 和其他工具来确定 CDR，但没有一个效果很好。因此，根据我们的数据库和实验结果，我们得出了自己的 CDR 定义。现在的挑战是根据定义将它们分开。

【解决方案2】：

我通过以下方法获得它，这对我来说非常适合找到 CDR1、2 和 3。我只需要定义 3 个不同的字典，其定义为前缀、后缀、最大引脚、固定位置，并将它们传递给以下代码。

在这里，我执行此操作以找到 CDR1，它为我提供了所需的输出。

dictionary_1={
            'cdr1_in': 'VCKASGFTFS',
            'cdr1_out':'WVRQAP',
            'max_pin_in':3,
            'position_c':2,
            'max_pin_out':2,
            'position_r':3
    
            }


def cdr_out(dictionary_1,x):
    count=0
    for i in range(len(x)-5):
        rider=x[i:i+6]
        # print(rider)
        if rider[2]=='R':
            # print(rider)
            for i in range(len(dictionary_1['cdr1_out'])):
                if rider[i]!=dictionary_1['cdr1_out'][i]:
                    count+=1
            if count<dictionary_1['max_pin_out']:
                # print(rider)
                return x.index(rider)

def cdr_in(dictionary_1,x):
    count=0            
    for i in range(len(x)-9):
        rider=x[i:i+10]
        # print(rider)
        if rider[1]=='C':
            # print(rider)
            for i in range(len(dictionary_1['cdr1_in'])):
                if rider[i]!=dictionary_1['cdr1_in'][i]:
                    count+=1
            if count<dictionary_1['max_pin_in']:
                # print(rider)
                y=x.index(rider)
                z=cdr_out(dictionary_1,x)
                cdr=x[y+10:z]
                return cdr
            
print(cdr_in(dictionary_1,x))

SYGMM

【讨论】：