具体实施方式

本发明提供了一种cDNA文库基因序列移码判断的方法，包括以下步骤：

1)将待比对的cDNA序列转换为Fasta格式、识别和去除接头序列，提取cDNA序列，得到输入序列；

2)构建包含候选蛋白的本地数据库；

3)使用blastx对所述步骤1)中的输入序列与所述步骤2)中的本地数据库进行比对，以匹配率最高的基因作为目标序列；

4)根据所述步骤1)中输入序列与所述步骤3)目标序列的位置比对信息，得到cDNA的移码情况；

所述cDNA的移码情况包括：

当所述输入序列比对到目标序列的起始位置为三的倍数加一时，为非移码；

当所述输入序列比对到目标序列的起始位置为三的倍数加二时，为移码一位；

当所述输入序列比对到目标序列的起始位置为三的倍数时，为移码两位；

步骤1)和步骤2)之间没有时间顺序限制。

本发明将待比对的cDNA序列转换为Fasta格式、识别和去除接头序列，提取cDNA序列，得到输入序列。

在本发明中，所述将待比对的cDNA序列转换为Fasta格式采用的软件优选的包括序列处理软件seqtk。

在本发明中，所述识别和去除接头序列采用的软件优选的包括awk的substr函数。

在本发明中，将待比对的cDNA序列转换为Fasta格式前，优选的还包括在将待比对cDNA序列转换为一行格式显示。在本发明中，将待比对cDNA序列转换为一行格式显示采用的软件包括序列处理软件seqtk。

本发明具体实施过程中，使用序列处理软件seqtk将所有待比对cDNA序列调整为每条序列一行显示，能够避免匹配接头的序列部分因换行导致无法匹配的情况，并将待处理的seq格式基因序列批量转化成Fasta格式。

本发明的方法可识别任意接头序列，并可批量提取接头后序列。具体步骤为使用awk的substr函数根据指定接头序列提取接头序列后第一个碱基之后(包含该第一个碱基)的所有序列(即只去掉接头序列及其前面部分的剩余序列，获得完整cDNA序列)，提取序列位置为Y＝X+L(X为接头序列第一个碱基的位置，L为接头序列长度)。

举例说明：

1)接头序列为ACAAGTTTGTACAAAAAAGTTGGX(SEQ ID NO.1，X为不固定碱基，可能是ATCG任意一个)，接头序列的长度为24(包含X)；

2)如果awk的substr函数获取的接头位置信息(即接头序列第一个碱基在整个序列上的位置)为80，则需要提取的序列起始位点为80+24＝104。整个序列的第104位为实际需要的CDS序列的起始碱基位置。

本发明构建包含候选蛋白的本地数据库。

在本发明中，优选的通过blast的makeblastdb进行本地数据库构建。在本发明中，所述候选蛋白的氨基酸序列根据比对需求进行选择。本发明对所述候选蛋白没有特殊限制，可以是某一物种、某些蛋白家族或者单一目标蛋白。本发明对所述本地数据库的数据来源没有特殊限制，来源于NCBI、uniprot及其他公开的蛋白数据库、个性化定制的蛋白序列或者校正及修正过的蛋白序列均可。本发明具体实施过程中，可构建个性化的本地数据库，进行精准比对。

得到本地数据库和输入序列后，本发明使用blastx对输入序列与本地数据库进行比对，输出的序列为目标序列。

在本发明具体实施过程中，在所研究物种为非常见物种，数据库中基因量少情况下，设置阈值为1e-5；在所研究物种为常见物种，基因组、转录组已被深度测序，数据库中录入的基因量多情况下，设置为1e-10。阈值越宽松，比对得到的结果越多，阈值越严格，比对得到的结果越少。低于上述阈值的序列不在输出结果中显示。

在本发明中，所述比对的参数包括：Identity、Gap、Align_length和E_value，根据各参数进行打分，打分最高的序列为输出的序列即目标序列。

得到目标序列后，本发明根据输入序列与所述目标序列的位置比对信息，判断cDNA的移码情况；

所述判断cDNA的移码情况包括：

当所述输入序列比对到目标序列的起始位置为三的倍数加一时，为非移码；

当所述输入序列比对到目标序列的起始位置为三的倍数加二时，为移码一位；

当所述输入序列比对到目标序列的起始位置为三的倍数时，为移码两位。

在本发明中，根据输入序列与所述目标序列的位置比对信息，优选的还包括确定输入序列和目标序列的匹配程度；所述确定输入序列和目标序列的匹配程度的判断标准包括：

比较输入序列的起始和终止位置、目标序列的起始和终止位置、gap和错配信息；

当输入序列的起始位置为1，终止位置为序列总长度；且目标序列的起始位置为1，终止位置为序列总长度，0错配，0gap，判定输入序列与目标序列完全匹配；

当输入序列的起始位置为1，终止位置小于序列总长度；且目标序列的起始位置为1，终止位置小于序列总长度，0错配，0gap，判定输入序列的5’端和目标序列完全匹配；

当输入序列的起始位置为1，终止位置为序列总长度；且目标序列的起始位置非1，终止位置为序列总长度，0错配，0gap，判定输入序列的3’端和目标序列完全匹配；

当输入序列的起始位置非1，终止位置小于序列总长度；且目标序列的起始位置非1，终止位置小于序列总长度，N错配，N gap，N大于等于0，且为整数，则输入序列与目标序列为不完全匹配，根据序列相似度算法输出具体匹配率。

下面将结合本发明中的实施例，对本发明中的技术方案进行清楚、完整地描述。显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

实施例1

按照下述步骤对cDNA序列进行比对：

(1)上传待比对序列及基因数据库到主机；

(2)将待比对序列格式转化为Fasta格式；

(3)将待比对序列与本地数据库中的参考基因或蛋白序列进行比对；

(4)根据候选基因与待比对基因匹配打分数值，输出与待比对基因相匹配的目的基因。比对结果如图1所示。

实施例2

按照下述步骤对cDNA序列进行比对，流程图如图2所示：

(1)上传待比对序列及基因数据库到主机；

(2)将待比对序列格式转化为Fasta格式；

(3)将待比对序列显示方式设置为1行；

(4)设置文库接头序列；

(5)在待比对序列中查找并删除接头上游及接头序列；

(6)将去除文库接头的cDNA序列，转换为氨基酸序列，在本地蛋白数据库中进行比对，分析匹配率；

(7)将比对打分最高的基因输出为最优结果，所述最优结果为目标基因；

(8)根据所设置的移码判断规则，计算待比对序列在载体中的表达是否移码；

(9)输出待比对cDNA在载体中是否移码及所比对到的目标基因详细信息。

对比实施例1

除了设置序列显示格式、设置文库接头序列、去除接头并计算待比对序列与目标基因第一匹配碱基位置数外与实施例2一致。

对比实施例2

增加序列显示格式，设置文库接头并去除接头序列，计算待比对序列与目标基因第一匹配碱基位置数，在实施例1快速批量比对的基础上，可以实现同时分析出cDNA在载体中是否移码。且对于不同的cDNA文库载体，仅需更换脚本中接头序列即可达到同样效果，简单快速。比对结果如图3所示，frameshift列所展示的即为是否移码的判断结果(黄色标亮所示)。其中frameshift列中的0表示cDNA在载体中没有移码，即插入基因正确表达；frameshift列中的1或2表示cDNA在载体中移码，需重新构建。

以上所述仅是本发明的优选实施方式，应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明原理的前提下，还可以做出若干改进和润饰，这些改进和润饰也应视为本发明的保护范围。