获取 vcf 文件中每个个体携带的突变数量
最近的工作涉及到处理大量种群数据，其中有个分析需要统计每个个体携带的突变数量，在网上搜寻了诸多方法后确定了一些比较高效的策略，固有此文记录。
本文的部分流程参考源：
BioStars "how to count variants par sample per chromosome in a vcf file?"
提取每个个体携带的突变数量该流程主要通过 plink2 实现：
12plink2 --vcf [vcf file] --out [plink output1]plink2 --pfile [plink output1] --sample-counts --out [plink output2]
[vcf file]：此处填入 vcf 文件的路径。
[plink output1]：此处填入 plink2 输出的文件前缀。
[plink output2]：此处填入 plink2 输出的文件前缀。
运行结束后，可以得到以 .scount 结尾的文件，其中各列的含义如下：
#IID：标识每个样本的编号。
HOM_REF_CT：纯合参考基因型计数，个体纯合参考基因型的位点数量 ...
vscode 在远程服务器上通过 jupyter notebook 使用 R
最近一段时间的关键分析要用到 R，但是办公电脑的机械硬盘前段时间出了问题，现在用 Rstudio 会非常卡顿，加上相关分析对计算资源的要求也比较高，就想着是否有什么办法能在服务器上也弄一个类似于 Rstudio 的可视化界面出来作为辅助。
不过看了一圈发现 Rstudio 的 server 版本似乎局限性不少，而且还要用到 root 权限。考虑到平时运行 Python 进行分析时大多数时候都在 jupyter notebook 上进行，就想着是否能套用过来，于是就有了这篇文章。
安装该文章的具体实现在 Linux + mamba 环境下完成，如果你没有 conda/mamba 环境，可以参考网上资料进行安装。
首先确认你的环境中有 R，如果没有可以通过以下命令安装：
12mamba install conda-forge::r-baseR # 确认可以进入 R 终端，使用 q() 退出
通过 mamba 安装 IRkernel 和 jupyter（使用 conda 也行，本文章以 mamba 为准）：
12mamba install r::r-irkernelmamba instal ...
簇突变的类型简介及识别方法（Clustered mutation classification）
基因组中，某些特殊的突变过程可能会导致一小段区域内连续的突变发生。这些突变在基因组中呈现成簇分布，因此被称作 clustered mutation（本文将称其为簇突变）。研究这些突变的分布及组成有助于我们揭示导致其发生的内源性和外源性过程。该文章将主要讲述簇突变的类别并列举一些会导致特定类别的生物学原因，同时介绍该领域中的某些生物信息学工具等。
簇突变类型簇突变类型主要被分作六类，其中五类针对碱基水平的变化，一类针对插入及缺失（InDel）：
DSB（Doublet Base Substitutions）：双碱基替换，表示这两个突变发生在相邻的碱基上。某些外源性过程例如紫外线损伤会导致双碱基替换的发生（CC>TT），此外 DNA 修复缺陷和聚合酶功能突变也会导致 DSB 发生。
MBS (Multiple Base Substitutions)：多碱基替换，表示在很短的序列范围内发生多个碱基突变且这些突变彼此相邻。由于该簇突变的出现数量很有限，因此尚未得到全面研究。
Omikli：源自希腊语，意为 雾 或 薄雾，也被称作 diffuse hypermutation（弥漫性超 ...
rsync 中排除文件的规则详解（rsync --exclude）
前段时间在备份服务器时，由于备份策略的变化，对 rsync 的 exclude patterns 进行了一些更加深入的探索。
这篇文章将介绍如何设计 rsync 的 exclude patterns 以排除特定目录或文件。
命令介绍rsync 是 Linux 和 Unix 中常用的文件复制和同步命令，相比 cp 命令而言更加高效和完善，其优点包括但不限于：
增量复制：rsync 只会复制那些自上次同步后更改的文件部分，极大地减少了数据传输量。例如 a 文件已经同步过一次并且没有发生变化，则 rsync 会跳过。同样，rsync 也会进行完整性检查以确定传输的数据无误。
保留属性：rsync 可以保持文件的权限、时间戳、软硬链接、所有权等属性。
灵活的过滤规则：支持使用包括和排除规则来选择同步的文件或目录。
排除规则rsync 通过参数 --exclude 或 --exclude-from 指定排除的文件模式。其中 --exclude-from 用于指定包含排除模式的文件（包含多个排除规则）。
排除指定文件输入文件名称（如 123.txt）即可排除所有命名为 123.txt 的文 ...
Snakemake pipeline 搭建的进阶教程
关于 SnakemakeSnakemake 是一个工作流管理系统，可用于创建易于迁移和复现的数据分析流程（pipeline），因此被广泛应用于生物学分析中。虽然 Snakemake 算是一个独立的编程语言，但其本质是基于 Python 搭建的，因此如果你对 Python 有一定了解，那么上手 Snakemake 就会更容易很多。并且 Python 技能的掌握对于某些情况下 Pipeline 下游分析部分的拓展是有必要的。
本文将结合多个近些年生物论文中公开的 Snakemake Pipeline，介绍如何搭建一个具有一定功能的 Pipeline。
Snakemake 安装：https://snakemake.readthedocs.io/en/stable/getting_started/installation.html
Snakemake Pipeline 示例教程：https://snakemake.readthedocs.io/en/stable/tutorial/basics.html
有利于提升阅读体验的一些条件：
对于 Python 编程已经具备一定的了解，包括文件 ...
对 bedgraph 文件进行 lowess 平滑操作
在进行一些特定的功能基因组学数据分析时，我们可能需要对 bedgraph 文件中每个 bin 的值进行一定的平滑操作，以降低随机噪声的影响并提供更好的可视化效果。例如：
Repli-seq / BrdU-seq 中量化得到的 Replication Timing
OK-seq / Pu-seq 中量化得到的 Replication Fork Directionality
以下是一个用于对 bedgraph 进行 lowess 平滑操作的 python script：
123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960import pandas as pdimport statsmodels.api as smimport argparsedef get_args(): parser = argparse.ArgumentParser(description='Perform lowess smooth on ...
获取基因组上 intron 区域及对应链信息的方法
流程前置条件：
在环境变量中可调用的 bedtools
本文参考：
Get intronic and intergenic sequences based on gff file from Biostars
https://www.biostars.org/p/112251/
本文可满足的需求：
得到特定区域的 bed 文件（例如 exon / intron 等）。
在得到区域信息的同时进行链信息的区分。
下载基因组注释文件（gtf）以人类最新版本的 gencode 注释为例：
1wget https://ftp.ebi.ac.uk/pub/databases/gencode/Gencode_human/release_46/gencode.v46.basic.annotation.gtf.gz
也可以选择已有的 gtf 文件进行。
得到 transcript 和 exon 区域信息12345678910111213141516GTF="/path/to/gencode.v46.basic.annotation.gtf.gz" # 改为自己实际的注释文件路径，建议使用压缩 ...
Ensembl VEP plugins 的使用方法（Alphamissense、dbNSFP 等）
前言以下是相关的官网链接，如果想要得到除本文内容以外更全面深入的了解，建议跳转相关链接进行查阅：
VEP plugins documentation:
https://grch37.ensembl.org/info/docs/tools/vep/script/vep_plugins.html
VEP install documentation:
https://grch37.ensembl.org/info/docs/tools/vep/script/vep_download.html
正文安装 VEP此处强烈建议使用 Docker 或者 Singularity 直接拉取 VEP 镜像进行分析，因为在 VEP image 中已经内置好所有的 plugin，无需自己手动下载。
考虑到 Docker 对权限的要求更高，因此这里以普适性更强的 Singularity 为例（Docker 用户有需求可通过文章开头链接前往官网进行参考）：
12# 拉取 VEP imagesingularity pull --name vep.sif docker://ensemblorg/ensembl- ...
fasterq-dump 下载 SRA 文件时报错的解决方法（err cmn_iter）
正文这几天下载 SRA，遇到的错误有：
fasterq-dump.3.1.0 err: cmn_iter.c cmn_read_uint8_array
fasterq-dump.3.1.0 err: cmn_iter.c cmn_read_String
以前使用 fasterq-dump 时只是偶然出现这些问题，重新下载也都能解决，但最近有些 SRA 一直下载失败，经过调查发现这种错误在大文件下载时出现异常频繁，且多次下载并没法有效解决问题，因此需要一个替代的方法。
解决方案Prefetch由于 fasterq-dump 直接通过 HTTP 下载得到 fastq 文件，该过程很可能由于一些问题中断从而导致下载失败。因此可以通过更稳定的 prefetch 先得到 sra 文件，再通过 fasterq-dump 提取 fastq 文件。
fasterq-dump fetches SRR on the fly via HTTP and there could be fatal errors during the transfer.prefetch eliminates transf ...
提取 bwa 比对中唯一比对 reads（uniquely mapped reads）的方法
本文由 Juse 基于以下资料进行撰写：
How to extract uniquely mapped reads from bam/sam produced by bwa-mem? from SEQanswers
Obtaining uniquely mapped reads from BWA mem alignment (filtering by q score does not seem to do the trick in my case) from Biostars
Obtaining uniquely mapped reads from BWA mem alignment from StackExchange
在此感谢 community 中各位大佬的无私分享。
正文在某些文章中，有时我们会看到作者提及到在后续分析中，ta 只考虑了 uniquely mapping reads（后文统称唯一比对）。从直觉上出发，使用唯一比对的思路不难理解 —— 这些比对更加精确，并使后续的分析具有更好的解释性。
目前来说，得到唯一比对的最主流方法是直接通过 MAPQ 进行过滤，但这并不 ...
bismark 分析 Bisulfite-Seq 数据的流程示例
前言以下是相关的官网链接和参考资料，如果想要得到除本文内容以外更全面深入的了解，建议直接跳转相关链接进行查阅：
Bismark github:
https://github.com/FelixKrueger/Bismark
Bismark documentation:
https://felixkrueger.github.io/Bismark/
本文比对及去重等流程均参考自 Bismark documentation
此外，以下内容在本文中不会涉及到，有需要请自行探索：
数据的质控
除 Bismark 以外其他软件的安装
通过以下流程可以得到：
基因组上 CpG 位点的甲基化信息
请阅读完全文后再进行相关分析，建议先使用 bismark 的示例数据观察相关命令是否存在问题。
正文下载最新版本的 bismark（至 2024/04/19）：
1$ wget https://github.com/FelixKrueger/Bismark/archive/refs/tags/v0.24.2.tar.gz
也可前往 github release 下载指定版本或可能已存在的更新 ...
通过 bowtie2 + macs3 分析 ChIP-Seq 数据
前言以下是相关的官网链接和参考资料，如果想要得到除本文内容以外更全面深入的了解，建议直接跳转相关链接进行查阅：
MACS3 github:
https://github.com/macs3-project/MACS
MACS3 documentation:
https://macs3-project.github.io/MACS/
本文比对及去重等流程参考自以下资料：
Alignment and filtering from Harvard Chan Bioinformatics Core (HBC)
Histone ChIP-seq pipeline from ENCODE project (Version 2.0 with bowtie2)
本文流程根据以下资料进行过优化：
sambamba -F “not duplicate” processed bam still have duplicated (sambamba issue#477) (Update 2024/4/24)
Encode blacklist (Update 2024/4/24)
以下内容在本文中不 ...
远程服务器无 root 权限时 LightGBM GPU 版的安装方法
正文首先感谢 LightGBM 团队的帮助：
https://github.com/microsoft/LightGBM/issues/6399
声明本文适用于在远程服务器（Linux）上安装 LightGBM 的情况。
不同的编译器版本在遇到报错时的解决方案可能略有不同，以下是我在编译时使用的相关版本：
12-- The C compiler identification is GNU 8.5.0-- The CXX compiler identification is GNU 4.8.5
本文默认已在服务器上配置有 CUDA，如果 CUDA 另有路径请通过 -DOpenCL_LIBRARY 和 -DOpenCL_INCLUDE_DIR 指定。
编译下载 LightGBM，并进行 CUDA 版本的编译：
123456git clone --recursive https://github.com/microsoft/LightGBMcd LightGBMmkdir buildcd buildcmake -DUSE_CUDA=1 .. # 指定 CUDA versionmake -j ...
全基因组 CpG (or nonCpG) 位点 bed 文件的获取方式
CpG 位点是基因组上的一种特殊序列，它在 DNA 的结构和功能中起着重要作用。其组成和甲基化影响着基因组上各区域的调控，是表观遗传学领域中重要的研究目标。
本文中的脚本主要用于提取基因组中 CpG 位点的位置，以方便一些着重于 CpG 位点的下游分析进行。此外在部分场景下，我们可能希望排除 CpG 位点进行分析（因为其独特的性质例如高突变率等），因此该脚本也可产生 AT&CG(nonCpG) 位点信息文件。
需要准备的文件和条件有：
基因组序列文件
Biopython 库（pip install biopython）
可能需要准备的前置条件：
已被放置于 $PATH 环境变量中的 bedtools
脚本内容：
1234567891011121314151617181920212223242526272829303132333435363738394041424344454647484950515253545556575859606162636465666768697071727374757677from Bio import SeqIOimport osimport ...