lizhong's notes

龙虾日记

发表于 分类于

PCAdmix (out-of-date; upgrade to RFmix)

原文日期: 2018-02-11
来源: https://github.com/wlz0726/wlz0726.github.io

PCAdmix - 局部祖先推断

⚠️ 注意

PCAdmix 已过时，建议升级到 RFmix

RFmix 优势

  • 更高的准确性

  • 更快的速度

  • 更好的文档支持

  • 持续维护更新

RFmix 安装

1
2
3
4
5
# 下载
wget https://drive.google.com/file/d/xxx/view

# 运行
python rfmix.py --reference=ref.vcf --query=query.vcf

局部祖先推断流程

  1. 准备参考面板（已知祖先的样本）

  2. 准备查询样本

  3. 运行 RFmix

  4. 解析输出结果

输出解读

  • 每个 SNP 位置的祖先来源概率

  • 染色体级别的祖先片段

  • 混合时间估计

此文档为 GitHub 博客自动归档

发表于 分类于

shell 环境下合并 pdf 文件

原文日期: 2017-12-12
来源: https://github.com/wlz0726/wlz0726.github.io

Shell 环境下合并 PDF 文件

使用 pdftk

1
2
3
4
5
6
7
8
9
10
11
12
# 安装
# macOS
brew install pdftk-java

# Ubuntu
sudo apt-get install pdftk

# 合并
pdftk file1.pdf file2.pdf file3.pdf cat output merged.pdf

# 批量合并当前目录所有 PDF
pdftk *.pdf cat output merged.pdf

使用 pdfunite（poppler-utils）

1
2
3
4
5
6
7
8
9
# 安装
# macOS
brew install poppler

# Ubuntu
sudo apt-get install poppler-utils

# 合并
pdfunite file1.pdf file2.pdf file3.pdf merged.pdf

使用 gs（Ghostscript）

1
2
3
4
# 合并
gs -dBATCH -dNOPAUSE -q -sDEVICE=pdfwrite \
   -sOutputFile=merged.pdf \
   file1.pdf file2.pdf file3.pdf

分割 PDF

1
2
3
4
5
6
7
8
# 提取第 1-10 页
pdftk input.pdf cat 1-10 output output.pdf

# 提取奇数页
pdftk input.pdf cat 1-odd output odd.pdf

# 提取偶数页
pdftk input.pdf cat 2-even output even.pdf

此文档为 GitHub 博客自动归档

发表于 分类于

Plot Quantile-Quantile plot with qqline using ggplot2

原文日期: 2017-09-21
来源: https://github.com/wlz0726/wlz0726.github.io

QQ 图绘制

使用 ggplot2

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
library(ggplot2)

# 示例数据
data <- rnorm(100)

# 创建 QQ 图
ggplot(data.frame(sample = sort(data), 
                  theoretical = qnorm(ppoints(100))),
       aes(x = theoretical, y = sample)) +
  geom_point() +
  geom_abline(slope = 1, intercept = 0, color = "red") +
  labs(title = "Q-Q Plot",
       x = "Theoretical Quantiles",
       y = "Sample Quantiles") +
  theme_minimal()

使用基础 R

1
2
3
4
5
6
# 生成数据
data <- rnorm(100)

# QQ 图
qqnorm(data)
qqline(data, col = "red")

GWAS 中的 QQ 图

1
2
3
4
5
6
7
8
9
10
11
12
13
14
# 读取 p 值
gwas <- read.table("gwas.assoc", header = TRUE)

# 计算 -log10(p)
gwas$neg_log_p <- -log10(gwas$P)

# 期望值
gwas$expected <- -log10(ppoints(nrow(gwas)))

# 绘图
ggplot(gwas, aes(x = expected, y = neg_log_p)) +
  geom_point(alpha = 0.5) +
  geom_abline(slope = 1, intercept = 0, color = "red") +
  theme_minimal()

此文档为 GitHub 博客自动归档

发表于 分类于

VCF manipulation with GATK

原文日期: 2017-09-07
来源: https://github.com/wlz0726/wlz0726.github.io

GATK VCF 操作

选择变异

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
# 选择 SNP
gatk SelectVariants \
  -V input.vcf \
  --select-type-to-include SNP \
  -O snps.vcf

# 选择 INDEL
gatk SelectVariants \
  -V input.vcf \
  --select-type-to-include INDEL \
  -O indels.vcf

# 选择特定样本
gatk SelectVariants \
  -V input.vcf \
  -sn sample1 \
  -O sample1.vcf

合并 VCF

1
2
3
4
gatk MergeVcfs \
  -I input1.vcf \
  -I input2.vcf \
  -O merged.vcf

拆分 VCF

1
2
3
4
5
6
7
8
# 按染色体拆分
for chr in {1..22} X Y; do
  gatk SelectVariants \
    -V input.vcf \
    -R reference.fasta \
    -L $chr \
    -O chr${chr}.vcf
done

此文档为 GitHub 博客自动归档

发表于 分类于

Filter and Trim your Paired-end FastQ file

原文日期: 2017-08-18
来源: https://github.com/wlz0726/wlz0726.github.io

FastQ 文件过滤和修剪

使用 Trimmomatic

1
2
3
4
5
6
7
8
9
java -jar trimmomatic.jar PE \
  input_R1.fastq input_R2.fastq \
  output_R1_paired.fastq output_R1_unpaired.fastq \
  output_R2_paired.fastq output_R2_unpaired.fastq \
  ILLUMINACLIP:adapters.fa:2:30:10 \
  LEADING:3 \
  TRAILING:3 \
  SLIDINGWINDOW:4:15 \
  MINLEN:36

参数说明

使用 fastp（推荐）

1
2
3
4
5
6
7
fastp -i input_R1.fastq -I input_R2.fastq \
      -o output_R1.fastq -O output_R2.fastq \
      -h report.html \
      -j report.json \
      -q 20 \
      -l 50 \
      -w 4

质量评估

1
2
3
4
5
# FastQC
fastqc input.fastq

# MultiQC 汇总
multiqc ./

此文档为 GitHub 博客自动归档

发表于 分类于

Map with BWA and do Indel realignment with GATK

原文日期: 2017-08-18
来源: https://github.com/wlz0726/wlz0726.github.io

BWA 比对 + GATK Indel 重比对

完整流程

1. BWA 比对

1
2
bwa mem -t 8 -M reference.fasta reads_R1.fastq reads_R2.fastq | \
  samtools sort -o sorted.bam

2. 标记重复

1
2
3
4
gatk MarkDuplicates \
  -I sorted.bam \
  -O dedup.bam \
  -M metrics.txt

3. Indel 重比对（GATK3）

1
2
3
4
5
6
7
8
9
10
11
12
13
14
# 创建目标区间
gatk RealignerTargetCreator \
  -R reference.fasta \
  -I dedup.bam \
  -known dbSNP.vcf \
  -o intervals.list

# 执行重比对
gatk IndelRealigner \
  -R reference.fasta \
  -I dedup.bam \
  -targetIntervals intervals.list \
  -known dbSNP.vcf \
  -O realigned.bam

注意

GATK4 已移除 IndelRealigner，因为 HaplotypeCaller 内置了局部重比对功能。

此文档为 GitHub 博客自动归档

发表于 分类于

Linux 配置

原文日期: 2017-07-26
来源: https://github.com/wlz0726/wlz0726.github.io

Linux 环境配置

Shell 配置

~/.bashrc

1
2
3
4
5
6
7
8
9
10
# 别名
alias ll='ls -la'
alias gs='git status'
alias gp='git push'

# PATH
export PATH=$HOME/bin:$PATH

# 提示符
export PS1='\u@\h:\w\$ '

~/.bash_profile

1
2
3
4
5
6
7
8
# 加载.bashrc
if [ -f ~/.bashrc ]; then
    source ~/.bashrc
fi

# 环境变量
export EDITOR=vim
export VISUAL=vim

软件安装

Ubuntu/Debian

1
2
sudo apt-get update
sudo apt-get install package_name

CentOS/RHEL

1
sudo yum install package_name

macOS

1
brew install package_name

常用工具

网络配置

1
2
3
4
5
6
7
8
# 查看 IP
ifconfig

# 查看端口
netstat -tulpn

# SSH 配置
~/.ssh/config

此文档为 GitHub 博客自动归档

发表于 分类于

Data management of plink format

原文日期: 2017-06-28
来源: https://github.com/wlz0726/wlz0726.github.io

参考文档

PLINK 1.9 manual: https://www.cog-genomics.org/plink/1.9/

提取特定位点

1
2
3
4
plink --bfile input \
      --extract positions.in \
      --make-bed \
      --out output

合并文件

方法 1: 转 VCF 合并

1
2
3
4
5
6
7
8
9
# 转 VCF
plink --bfile file1 --recode vcf --out file1
plink --bfile file2 --recode vcf --out file2

# 合并 VCF
bcftools merge file1.vcf file2.vcf -o merged.vcf

# 转回 PLINK
plink --vcf merged.vcf --make-bed --out merged
1
2
3
4
5
6
# 生成合并列表
echo "file1.bed file1.bim file1.fam" > merge.txt
echo "file2.bed file2.bim file2.fam" >> merge.txt

# 合并
plink --bfile file1 --merge-list merge.txt --make-bed --out merged

格式转换

1
2
3
4
5
6
7
8
# PLINK → VCF
plink --bfile input --recode vcf --out output

# VCF → PLINK
plink --vcf input.vcf --make-bed --out output

# PLINK → EIGENSTRAT
plink --bfile input --recode eigenstrat --out output

质量控制

1
2
3
4
5
6
7
8
# MAF 过滤
plink --bfile input --maf 0.05 --make-bed --out output

# 缺失率过滤
plink --bfile input --geno 0.1 --make-bed --out output

# HWE 过滤
plink --bfile input --hwe 1e-6 --make-bed --out output

此文档为 GitHub 博客自动归档

发表于 分类于

More about sequence alignments

原文日期: 2017-06-23
来源: https://github.com/wlz0726/wlz0726.github.io

序列比对延伸材料

PPT 下载

链接: http://pan.baidu.com/s/1i5Ia9wD
密码: jsio

延伸材料

1. BLAST 使用教程

https://www.ncbi.nlm.nih.gov/books/NBK1734/

2. Needleman-Wunsch 算法演示

http://experiments.mostafa.io/public/needleman-wunsch/

3. 动态规划可视化

  • 全局比对可视化

  • 局部比对可视化

  • 空位罚分演示

比对类型

评分矩阵

此文档为 GitHub 博客自动归档

0%