上一期我们对新冠变异株核酸检测的方法原理进行了简介,本期我们主要讲一下,如何利用python语言实现我们特异性检测的目标。

一、文件读取

接着上一期我们下载新冠Alpha、Beta、Delta、Gamma株的全部基因组序列后,我们发现同一株型的新冠基因组,以“>”分开不同测序基因组,序列总共有上千条,为了后面程序简便,我们将同一株型的新冠基因组文件读取并储存数组中,如Alpha株,其数组形式为:Alpha[“序列1”,“序列2”,“序列3”,......,“序列n”],其余毒株方法相同,代码如下:

Alpha = []
i = 0
a = ""
f = open("/home/lxh/Alpha.fasta","r")
for line in f.readlines():if ">" in line:i += 1Alpha.append(a)a = ""if ">" not in line:a += line.replace("\n","")del Alpha[0]f.close()

简单说明下,上面代码意思是逐行读取,当遇到“>”建立新的数组元素,没有遇到的时候也就是将每一行核酸序列末尾的换行符“\n”去掉,加在一起变成一个字符串,这样一个字符串成为毒株数组中的一个元素,以此类推,完成所有基因组的读取和写入。

二、引物设计

本次引物长度暂定为20bp,如要设计Alpha毒株的引物,那么引物必须满足能够同时出现Alpha数组的所有元素中,但不出现在其他毒株元素中。我们参考宏基因组测序kermer的分析方法。假设有30个序列,如果引物长度为20bp,那么会出现30-20+1=11条连续而且大概率不相同的引物。代码如下:

Alpha_kermer = []
for m in range(i-1):n = len(Alpha[m])for mm in range(n-19):Alpha_kermer.append(Alpha[m][mm:(20 + mm)])

注意循环语句最后一句一定要从头空一行,表示循环语句从这里结束,新手一定注意,不懂的可以搜索WC3School,里面有各组编程语言的语法教程。

以上代码将Alpha所有的基因组裁剪并形成20bp的数组,因为这些元素来自于同一株的新冠基因组,所以理论上就会出现大量相同序列的元素,下面的代码将数组中元素进行了频次统计,并进行元素出现频次高低的排序。新形成的bb成为一个无重复,频率由高到低出现的数组,通过检验bb的第一个元素bb[0]为一串polyA,果断放弃,第二个元素bb[1]为正常序列,出现频率为3879次,而我们上面写入的新冠基因组也有3879次,这说明此序列出现在每个基因组中。代码如下:

from collections import Counter
collection_words = Counter(Alpha_kermer)
ibb = sorted(collection_words.items(),key = lambda x:x[1],reverse = True)

出现3879次的序列远不止这一条,我们将出现频次为3879次的所有序列检索并形成新的数组,为了减少意外情况出现和重复计算,我们将结果导出保存,以便下次直接读取应用,代码如下:

f = open("Alpha_kermer.txt","a+")
fk = open("Alpha_kermer_all_in.txt","a+")
k_number = 0
for key,value in bb:f.write(key + "\n")if value == (i - 1):k_number += 1fk.write(key + "\n")f.close()
fk.close()
k_number

文件"Alpha_kermer.txt"为出现在Alpha中长度为20bp的所有引物序列,文件Alpha_kermer_all_in.txt为同时出现在所有Alpha株中长度为20bp的引物序列,两者为包含关系:“或”、“和”关系。

其他株的代码,只需要查找和替代Alpha,换成其他株名就可以一键copy and run。

新冠病毒变异株核酸检测引物设计——代码实现1相关推荐

  1. 港科夜闻|香港科大新研究显示人体T细胞免疫反应可有效应对新冠病毒变异株Omicron...

    关注并星标 每周阅读港科夜闻 建立新视野 开启新思维 1.香港科大新研究显示人体T细胞免疫反应可有效应对新冠病毒变异株Omicron.透过接种新冠疫苗或经感染新冠肺炎后所产生的T细胞,能有效清除被感染 ...

  2. 实验一 词法分析程序设计_重庆新增一生物2级安全实验室,将着眼新冠病毒等药物检测分析...

    5月12日,上游新闻·重庆商报记者获悉,位于两江新区的迪纳利医药科技有限责任公司(下简称"迪纳利")将建立生物2级安全实验室(简称"BSL-2"),该实验室主要 ...

  3. 阿联酋研发新冠病毒快速激光检测技术

    阿联酋阿布扎比--(美国商业资讯)--在阿布扎比证券交易所上市的International Holdings Company (IHC)旗下药物研究分支机构QuantLase Imaging Lab宣 ...

  4. 关于一种新的空气内新冠病毒检测方式的诸多设想

    关于一种新的空气内新冠病毒检测方式的诸多设想 引言 目前空气中新冠病毒的检测方式大多为首先利用特殊吸附膜将气溶胶态的病毒富集,之后溶解为液态以进行常规的核酸检测1,而该方法由于涉及到的步骤较多,操作复 ...

  5. 最快60秒完成新冠病毒核酸对比 阿里云向社会免费开放基因计算服务

    全球疫情肆虐,各大科技公司都在竭尽全力抗击疫情.3月13日,阿里云对外宣布,将向医疗科研机构.疾控中心等一线病毒研究机构免费开放基因计算服务,可大幅提升宏基因组测序.疫苗研发相关的处理效率,最快只需6 ...

  6. 一文读懂测序技术在新冠病毒检测中的应用(文末附FAQ)

    来源:生物探索 随着世界疫情的发展,多个国家进入公共卫生紧急状态,全球科学家都在抓紧研究更好的检测.治疗.防控手段.从最初未知β属冠状病毒的快速鉴定到病毒序列的完整破译,再到病毒序列的变异监测,高通量 ...

  7. 新风医疗集团就私有化交易达成最终合并协议;​康泰生物成功分离新冠德尔塔变异株 | 医药健闻...

    | 行业焦点 新风医疗集团就私有化交易达成最终合并协议.新风医疗集团宣布已就私有化与 Unicorn II Holdings Limited及其下属公司签署了最终合并协议.本次交易公司的股权价值约为 ...

  8. qPCR与新冠病毒核酸检测

    qPCR与新冠病毒核酸检测 这篇文章简单记录一下qPCR核酸定量的原理以及它在新型冠状病毒检测中的应用,分享给身边好奇核酸检测的朋友,参考视频链接:https://www.bilibili.com/v ...

  9. 京东健康上线“新冠病毒”核酸检测 在线预约服务

    近日,京东健康与北京金域医学实验室达成合作,上线"新冠病毒"核酸检测的在线预约服务,成为全国首个提供核酸检测服务在线下单及预约的平台. 据介绍,用户通过京东APP搜索"核 ...

最新文章

  1. 修改密码后服务器断开连接,SSH无需密码登录服务器且保持连接不断开的方法
  2. JS设置cookie、读取cookie、删除cookie
  3. 文章章节序号编排常识
  4. webpack4.x最详细入门讲解
  5. 文献记录(part69)--公平性机器学习中基于分类间隔的歧视样本发现和消除算法
  6. P3768 简单的数学题 [狄利克雷卷积,杜教筛,莫比乌斯反演]
  7. 任务并行VS数据并行
  8. LeetCode 106. 已知中序后序 求二叉树
  9. 建议你一定要尝试的副业排名TOP1
  10. js页面传值php页面,php实现跳转传值有什么方法,js页面跳转传值
  11. 开源地图引擎openlayers_由quot;地图quot;到quot;指南针quot;:疫后智能营销的演化逻辑...
  12. java系统性能优化之mysql数据库优化
  13. Linux源码安装pgadmin4,赵彦昌博客 - linux ubuntu 安装pgadmin4
  14. B2C电子商务开发的网店管理系统
  15. 【生物信息学】——Metagenomics宏基因组学分析流程浅谈
  16. 罗技鼠标的蓝牙适配器无效的解决办法
  17. what to benefit from the C++14 Standard
  18. video.js视频播放器进度条标记的功能实现
  19. 双休和单休区别大吗?
  20. Huawei 5G MiFi E6878-370 VS E6878-870

热门文章

  1. cadence allegro导入dxf文件
  2. Qt --- QTreeWidget 树形控件实例遇到的问题
  3. Java 基础巩固:内部类的字节码学习和实战使用场景
  4. UltraEdit打开文件,中文显示为乱码的解决方法
  5. LaTeX 制作(跨页)长表格
  6. 发泡聚苯乙烯(EPS)的全球与中国市场2022-2028年:技术、参与者、趋势、市场规模及占有率研究报告
  7. 2008-2020年数据上市公司高管团队异质性数据包含Stata代码
  8. 过拟合解决方案 —— early stopping
  9. html 加减法,加减法速算技巧太神奇了!
  10. 浅谈运营商网络业务限速(上)