windows ubuntu 子系统:肿瘤全外篇,2. fq 数据质控,比对。

首先我们先下载一组全外显子测序数据。nabi sra库,随机找了一个。

来自受试者“16177_CCPM_1300019”(SRR28391647, SRR28398576)的样本“16177_CCPM_1300019_BB5”的基因组DNA配对端文库“0369547849_Illumina_P5-Popal_P7-Hefel”的Illumina随机外显子测序

下载下来,转为两个配对的fq文件。过程可参考3.windows下Ubuntu,sratoolkit软件,从ncbi的sra数据库下载数据。_sratools下载数据-CSDN博客

这样我们得到了两个配对的fq文件,如果太大,可以压缩一下。

1.质控

我们直接进行fastp进行质量控制,这时我们可以实现以下流程化。

 ls *_1.fq.gz >1
 ls *_2.fq.gz >2
paste 1 2 >  merge
sed -i "s/.fq.gz//g" merge

cat merge | while read id ; do
    fastp -i "${id}_1.fq.gz" -I "${id}_2.fq.gz" -o fastp/"${id}_clean_1.fq.gz" -O fastp/"${id}_clean_2.fq.gz" -j fastp/"${id}.json" -h fastp/"${id}.html"
done

这样我们就可以对一个目录下的样本进行批量指控了,将它们填写到脚本中,可以直接运行。

2.比对并排序

接下来,我们使用bwa mem模式对经过质控的fq文件进行比对,生成bam文件,进行排序。

cat ../merge | while read id ; do
    bwa mem -t 20 -R "@RG\tID:${id}\tLB:${id}\tPL:Illumina\tSM:${id}" /mnt/h/db/bwa.db/hg38.fa ./"${id}_clean_1.fq.gz" ./"${id}_clean_2.fq.gz" | samtools sort -@ 2 -o human/"${id}.bam"
done

        这上面的-R参数很重要。-R 接的是 Read Group的字符串信息,这是一个非常重要的信息,以@RG开头,它是用来将比对的read进行分组的。不同的组之间测序过程被认为是相互独立的,这个信息对于我们后续对比对数据进行错误率分析和Mark duplicate时非常重要。在Read Group中,有如下几个信息非常重要:

(1) ID,这是Read Group的分组ID,一般设置为测序的lane ID(不同lane之间的测序过程认为是独立的),下机数据中我们都能看到这个信息的,一般都是包含在fastq的文件名中;

(2) PL,指的是所用的测序平台,这个信息不要随便写!特别是当我们需要使用GATK进行后续分析的时候,更是如此!这是一个很多新手都容易忽视的一个地方,在GATK中,PL只允许被设置为:ILLUMINA,SLX,SOLEXA,SOLID,454,LS454,COMPLETE,PACBIO,IONTORRENT,CAPILLARY,HELICOS或UNKNOWN这几个信息。基本上就是目前市场上存在着的测序平台,当然,如果实在不知道,那么必须设置为UNKNOWN,名字方面不区分大小写。如果你在分析的时候这里没设置正确,那么在后续使用GATK过程中可能会碰到类似如下的错误:

 ERROR MESSAGE: The platform (xx) associated with read group GATKSAMReadGroupRecord @RG:xx is not a recognized platform.

这个时候你需要对比对文件的header信息进行重写,就会稍微比较麻烦。

我们上面的例子用的是PL:illumina。如果你的数据是CG测序的那么记得不要写成CG!而要写COMPLETE。

(3) SM,样本ID,同样非常重要,有时候我们测序的数据比较多的时候,那么可能会分成多个不同的lane分布测出来,这个时候SM名字就是可以用于区分这些样本。

(4) LB,测序文库的名字,这个重要性稍微低一些,主要也是为了协助区分不同的group而存在。文库名字一般可以在下机的fq文件名中找到,如果上面的lane ID足够用于区分的话,也可以不用设置LB;

除了以上这四个之外,还可以自定义添加其他的信息,不过如无特殊的需要,对于序列比对而言,这4个就足够了。这些信息设置好之后,在RG字符串中要用制表符(\t)将它们分开。

3.标记PCR重复,使用picard

cat ../merge | while read id ; do
    java -jar /mnt/h/softwore/picard/picard.jar MarkDuplicates \
    I="${id}.bam" \
    O="${id}.markup.bam" \
    M="${id}.markdup_metrics.txt"
done

  1. cat ../merge: 这个命令会将 merge 文件的内容输出到标准输出流。
  2. while read id; do ... done: 这是一个 while 循环,它会逐行读取 cat 命令的输出,并将每一行的内容赋值给变量 id
  3. java -jar /mnt/h/softwore/picard/picard.jar MarkDuplicates ...: 这是 Picard 工具的命令,用于标记重复的 reads。具体的参数解释如下:
    • -jar /mnt/h/softwore/picard/picard.jar: 指定 Picard 工具的 JAR 文件路径。
    • MarkDuplicates: 执行标记重复 reads 的操作。
    • I="${id}.bam": 输入 BAM 文件的路径和文件名,${id}.bam 表示根据当前循环的 id 变量构建的 BAM 文件名。
    • O="${id}.markup.bam": 输出标记重复 reads 后的 BAM 文件的路径和文件名,${id}.markup.bam 表示根据当前循环的 id 变量构建的标记后的 BAM 文件名。
    • M="${id}.markdup_metrics.txt": 输出标记重复 reads 后的统计信息文件的路径和文件名,${id}.markdup_metrics.txt 表示根据当前循环的 id 变量构建的统计信息文件名。

4.samtools建立索引

cat ../merge | while read id ; do
    samtools index "${id}.markup.bam" -@ 10
 done

接下来,我们就可以进入GATK4的流程了。

        这是我最近跑的一个流程,因为我在很久之前就跑过,所以这次没有很多基础信息,我稍后会分享一个资料,来对这些流程进行一个原理的解释。这个流程就是跑一下主要过程,想要建立一个流程,得要考虑很多东西,选择软件,修改添加参数啥的等,还是得有具体的项目实施。

本文来自互联网用户投稿,该文观点仅代表作者本人,不代表本站立场。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如若转载,请注明出处:/a/552369.html

如若内容造成侵权/违法违规/事实不符,请联系我们进行投诉反馈qq邮箱809451989@qq.com,一经查实,立即删除!

相关文章

SGI_STL空间配置器源码剖析(一)总览

SGI 全称为 Silicon Graphics [Computer System] Inc. 硅图[计算机系统] 公司,SGI_STL是SGI实现的C的标准模板库。 SGI STL的空间配置器包括一级和二级两种。 一级空间配置器allocator采用malloc和free来管理内存,这与C标准库中提供的allocator是相似的…

VS集成vcpkg

VS集成vcpkg 下载vcpkg 下载vcpkg git clone https://github.com/Microsoft/vcpkg.git安装vcpgk,文件目录 .\bootstrap-vcpkg.bat集成到vs2022中 # 集成到项目 vcpkg integrate project vcpkg integrate installPS C:\Users\Administrator> vcpkg integrate…

大模型开发轻松入门——(1)从搭建自己的环境开始

pip install openai import openai import osfrom dotenv import load_dotenv, find_dotenv _ load_dotenv(find_dotenv())openai.api_key os.getenv(OPENAI_API_KEY)

如何选择投资交易策略?很简单,只需回答fpmarkets6个问题

刚迈出交易的第一步的投资新手们,是不是还没有选择策略?外汇市场上的交易策略是一种算法,可以让投资者以最低的风险尽快实现目标。目标通常是获得一定比例的利润。 那么如何选择投资交易策略?很简单,只需回答fpmarkets…

计算机网络 2.2数据传输方式

第二节 数据传输方式 一、数据通信系统模型 添加图片注释,不超过 140 字(可选) 1.数据终端设备(DTE) 作用:用于处理用户数据的设备,是数据通信系统的信源和信宿。 设备:便携计算机…

酒店餐厅装水离子雾化壁炉前和装后对比

酒店餐厅装水离子雾化壁炉前和装后的对比可以体现出餐厅氛围和客户体验的显著改变: 装前: 普通的氛围:餐厅可能显得比较普通,缺乏特色或独特的装饰元素。 视觉上缺乏焦点:餐厅空间可能显得相对平淡,缺乏…

如何在MacOS上使用OpenHarmony SDK交叉编译?

本文以cJSON三方库为例介绍如何通过OpenHarmony的SDK在Mac平台进行交叉编译。 环境准备 SDK准备 我们可以通过 openHarmony SDK 官方发布渠道下载对应mac版本的SDK,当前OpenHarmony MAC版本的SDK有2种,一种是x86架构,另一种是arm64&#x…

【小白学机器学习13】一文理解假设检验的反证法,H0如何设计的,什么时候用左侧检验和右侧检验,等各种关于假设检验的基础知识

目录 前言: 目标 1 什么叫 假设检验 1.1 假设检验的定义 1.1.1 来自百度百科 1.1.2 维基百科 1.2 假设检验的最底层逻辑:是反证法思想 1.3 假设检验的底层构造:小概率反证法思想 2 什么叫反证法 2.1 反证法的概念 2.1.1 来自百度…

HarmonyOS开发实例:【任务延时调度】

介绍 本示例使用[ohos.WorkSchedulerExtensionAbility] 、[ohos.net.http]、[ohos.notification] 、[ohos.bundle]、[ohos.fileio] 等接口,实现了设置后台任务、下载更新包 、保存更新包、发送通知 、安装更新包实现升级的功能。 效果预览 使用说明 安装本应用之…

基于Python的机器学习的文本分类系统

博主介绍:✌程序员徐师兄、7年大厂程序员经历。全网粉丝12w、csdn博客专家、掘金/华为云/阿里云/InfoQ等平台优质作者、专注于Java技术领域和毕业项目实战✌ 🍅文末获取源码联系🍅 👇🏻 精彩专栏推荐订阅👇…

PyTorch深度学习入门-2

PyTorch深度学习快速入门教程(绝对通俗易懂!)【小土堆】_哔哩哔哩_bilibili 一、神经网络的基本骨架 --nn.Module Neutral network torch.nn — PyTorch 2.2 documentation * import torch from torch import nnclass xiaofan(nn.Module):…

探索未来:人工智能—图像分类的发展与核心技术

引言 在当今数字化时代,图像已经成为我们生活中不可或缺的一部分,而人工智能技术的发展为图像处理和分析提供了巨大的机遇和挑战。其中,图像分类作为人工智能领域的一个重要应用,在诸多领域中发挥着关键作用。 人工智能在图像分类…

Pascal VOC(VOC 2012、VOC 2007) 数据集的简介

一、数据集介绍 PascalVOC(2005~2012)数据集是PASCAL VOC挑战官方使用的数据集。该数据集包含20类的物体。每张图片都有标注,标注的物体包括人、动物(如猫、狗、岛等)、交通工具(如车、船飞机等)、家具(如椅…

多线程意义

直接上代码 我们来看两个程序 由一个线程和两个线程运行的区别&#xff1a; 单线程&#xff08;main&#xff09;&#xff1a; public static void test(){long a 0;long b 0;for(long i 0; i < 10000000000l; i){a;}for(long i 0; i < 10000000000l; i){b;}} 多…

MySQL Prepared语句(Prepared Statements)

在数据库应用中&#xff0c;很多SQL语句都会重复执行很多次&#xff0c;每次执行可能只是where条件中的变量值不同&#xff0c;但MySQL依然会解析SQL语法并生成执行计划。对于这类情况&#xff0c;可以利用prepared语句来避免重复解析SQL的开销。 文章目录 一、prepared语句优…

蓝桥杯(基础题)

试题 C: 好数 时间限制 : 1.0s 内存限制: 256.0MB 本题总分&#xff1a;10 分 【问题描述】 一个整数如果按从低位到高位的顺序&#xff0c;奇数位&#xff08;个位、百位、万位 &#xff09;上 的数字是奇数&#xff0c;偶数位&#xff08;十位、千位、十万位 &…

《系统分析与设计》实验-----在线书店系统 需求规格说明书 哈尔滨理工大学PLUS完善版

文章目录 需求规格说明书1&#xff0e;引言1.1编写目的1.2项目背景1.3定义1.4参考资料 2&#xff0e;任务概述2.1目标2.2运行环境2.3条件与限制 3&#xff0e;数据描述3.1静态数据3.2动态数据3.3数据库介绍3.4数据词典3.5数据采集 4&#xff0e;功能需求4.1功能划分4.2功能描述…

ES-全文搜索

模糊查询&#xff1a; 写数据通过id路由到master分片 查询数据到一个节点&#xff0c;该节点会作为一个调度节点判断负载等情况将请求转发到真正节点&#xff08;一般し轮询&#xff09;

C语言-指针

1. 指针是什么 指针理解的2个要点&#xff1a; 1.1. 指针是内存中一个最小单元的编号&#xff0c;也就是地址 1.2 平时口语中说的指针&#xff0c;通常指的是指针变量&#xff0c;是用来存放内存地址的变量 总结&#xff1a;指针就是地址&#xff0c;口…

vue+element作用域插槽

作用域插槽的样式由父组件决定&#xff0c;内容却由子组件控制。 在el-table使用作用域插槽 <el-table><el-table-column slot-scope" { row, column, $index }"></el-table-column> </el-table>在el-tree使用作用域插槽 <el-tree>…