●在fastq-dump拆分SRA文件時遇到報錯

image.png

解決方案：

因為NCBI上的下載鏈接從http變?yōu)榱薶ttps，所以安裝最新版sratoolkit即可解決問題

●使用fasterq-dump拆分SRA文件,速度更快

fasterq-dump --split-files SRR934398.sra -e 10 -p

##--split-files最好是--split-3，但是有些特殊情況--split-3拆分不開雙端數(shù)據(jù)就用--split-files。

●遇到報錯error while loading shared libraries: libcrypto.so.1.0.0: cannot open shared object file: No such file or directory

1、先查看libcrypto.so.1.0所在目錄

image.png

發(fā)現(xiàn)沒有libcrypto.so.1.0.0這個文件，因此，嘗試建立一個叫libcrypto.so.1.0.0的符號鏈接向/usr/lib64/libcrypto.so.1.0.2k

2、建立軟鏈接

ln -s /usr/lib64/libcrypto.so.1.0.2k ~/miniconda3/lib/libcrypto.so.1.0.0

●統(tǒng)計測序深度，使用vcf文件

vcftools --vcf test.vcf \
         --depth -c \
         > depth_summary.txt

●提取染色體片段

       vcftools --gzvcf Duroc.vcf.gz \
                --chr NC_010484.4 \
                --from-bp 181225 \
                --to-bp 182187 \
                --out MC1R-Duroc.vcf \
                --recode \
                --recode-INFO-all 

●提取文件中的某幾列

awk '{print$1,$2}' fileame.vcf > filename.txt
##$1代表第一列

●根據(jù)位置提取vcf文件對應(yīng)位點(diǎn)的信息

vcftools --vcf 12_28.filter.snps.indels.vcf --positions test.txt --out test --recode

●提取某一列數(shù)值滿足條件的列

awk -F'\t' '{if ($3 == 1) print $1\t$2\t$3}' t.txt > 1.txt

#以Tab鍵分割
awk 'BEGIN{IFS='\t'}{if ($5 > 0) print $1,$2,$2+1,$5}' 12.26-CHINA.freq.frq > 12.26-CHINA.maf0.frq.txt

●提取某些樣本

bcftools view -S id.txt 20211005_sheep_222_total.vcf.gz > tibetan_36.vcf 

###●其中 id.txt 為一列樣本id

●去除vcf文件中帶*的等位基因

grep -v "*" JBC-geno005-maf005.vcf > JBC-geno005-maf005.filter.vcf

●格式轉(zhuǎn)化

1 bed、bim、fam轉(zhuǎn)vcf

plink --allow-extra-chr \
      --chr-set 26 \
      -bfile xll \
      --recode vcf-iid \
      --out xll

2 bed、bim、fam轉(zhuǎn)map、ped

plink --allow-extra-chr \
      --chr-set 26 \
      -bfile filename \
      --recode \
      --out filename

3 ped、map轉(zhuǎn)bed、bim、fam

plink --allow-extra-chr \
      --chr-set 26 \
      --file tibetan_36 \
      --make-bed \
      --out tibetan_36

4 map、ped轉(zhuǎn)為vcf

plink --allow-extra-chr \
      --chr-set 26 \
      -file xll \
      --recode vcf-iid \
      --out xll

5 vcf轉(zhuǎn)ped、map

plink --allow-extra-chr \
      --chr-set 26 \
      --vcf tibetan_36.vcf \
      --recode \
      --double-id \   
      --out tibetan_36
###double-id兩個family id 和idividual id一樣，所以加上了這個代碼

6 vcf轉(zhuǎn)bed、bim、fam

plink --allow-extra-chr \
      --chr-set 26 \ 
      -vcf XXX.vcf \
      --make-bed \
      --double-id \
      --out XXX

●缺失率統(tǒng)計

## 按照位點(diǎn)統(tǒng)計
vcftools --gzvcf test.vcf.gz \
         --missing-site \
         --out test.SNP_missing 
## 按照個體統(tǒng)計
vcftools --vcf  test.vcf   \
         --missing-indv \
         --out test.SNP_missing

●過濾indel和snp

## INDEL
vcftools --remove-indels \
         --recode \
         --recode-INFO-all \
         --vcf test.vcf \
         --stdout \
         > test.snp.vcf
## SNP
vcftools --keep-only-indels  \
         --recode \
         --recode-INFO-all \
         --vcf test.vcf \
         --stdout \
         > test.indel.vcf

●去除多等位基因及indel

bcftools view -m 2 \
              -M 2 \
              --type "snps"  test.vcf.gz \
              -Ov \
              -o test.record.snps.vcf.gz

## 注意一下：-O為輸出文件的格式，其中z為壓縮的vcf文件，v為正常的vcf文件，
vcftools --vcf   test.vcf \
         --remove-indels \
         --min-alleles 2 \
         --max-alleles 2 \
         --recode \
         --recode-INFO-all \
         --stdout \
         > test.miss.snp.vcf