plink commands

convert-formats.sh

# ------------------
# from bgen to bed
# ------------------
plink2 --bgen myfile.bgen \
  # filter by chr or by snps names
  --chr 1 --extract {params.snps} \ 
  --make-bed --out myoutput
# 1. plink will create temporary/intermediate files such as *-temporary.{pgen,psam,pvar}
# 2. sample ids in .fam will be `ID1_ID1`

#---------------------
# from bed to vcf
#---------------------
plink --bfile $bfile --chr $chr  \
  --keep-allele-order --recode vcf-iid --out $chr \
  --threads $cores

bgzip -f $chr.vcf

tabix -f $chr.vcf.gz

prepare-data-phasing.sh

# Tractor pipeline: https://github.com/eatkinson/Post-QC/blob/master/CohortDataQC_final.sh
# Plink command example to merge two bed files: https://www.biostars.org/p/101191/#101247

# Annotations.
# rsid names from snpdb
# https://hgdownload.cse.ucsc.edu/goldenPath/hg38/database/
GENOME=hg38
SNPBUILD=144
curl -s http://hgdownload.cse.ucsc.edu/goldenPath/$GENOME/database/snp$SNPBUILD.txt.gz | zcat | cut -f 2,3,4,5,6,7,10,16 > dbsnp.$SNPBUILD.$GENOME.bed

qc-plink.sh

#---------------------
# remove duplicates
#---------------------
# - https://www.cog-genomics.org/plink/1.9/data#list_duplicate_vars

plink --bfile $bfile --chr $chr \
  --list-duplicate-vars suppress-first ids-only
# - output: plink.dupvar
plink --bfile $bfile --chr $chr  \
  --exclude plink.dupvar \
  --make-bed --out out
# - output: out.bed

run-bcftools.sh

# References
# - https://gist.github.com/elowy01/93922762e131d7abd3c7e8e166a74a0b

#---------
# stats
#---------
# number of samples
bcftools query -l $chr.vcf | wc -l
# number of variants
bcftools query -f '%POS\n' $chr.vcfist | wc -l

#----------
# filter
#----------
# filter by variant IDs
# - https://github.com/samtools/bcftools/issues/776
bcftools view -i'ID=@test.snp' my.vcf 
# filter by sample IDs
bcftools view -s sample1,sample2 file.vcf > filtered.vcf
bcftools view -S sample_file.txt file.vcf > filtered.vcf

#----------------
# update fields
#----------------
# https://www.biostars.org/p/363263/
# replace variant ID by a custom name 'chr:pos:ref:alt'
bcftools annotate --set-id '%CHROM:%POS:%REF:%FIRST_ALT' input.vcf > output.vcf

#----------
# index
#----------
bgzip -f $chr.vcf
tabix -f $chr.vcf.gz

run-plink.sh

# -----------------------------
# merge bed files across chr.
# -----------------------------
for i in {1..22}; do
  echo chr$i >> mergelist.txt # bfile, e.g. without .bed extension
done
plink --merge-list mergelist.txt --make-bed --out auto # resulted file: auto.bed