###### tags: `解析` `NGS` `ChIP-seq` `omni-ATAC` `EM-seq` `shell script`

NGS Data Analysis (QC 〜 Mapping)
===

　いちおう完成版？( 2019/12/6 バグ等 何か問題がある場合はマエダまで...)

### 0. インストール
☆ Skewer、PRINSEQ、FastQC が必要
<br>

**0-1-1. Skewerのダウンロード**
*ダウンロードしたいディレクトリに移動してダウンロード*
　`git clone https://github.com/relipmoc/skewer.git `

*skewerディレクトリに移動*
　`cd skewer/`

　`make`
<br>

**0-1-2. SkewerのPATHを通す**
`.bash_profile` に vi,vim,echo等で インストールしたskewerの場所を加え、変更を反映
*（例）*
　`PATH=~/ProgramFiles/skewer:$PATH`
*変更したbashprofileを反映*
　`source ~/.bashrc`
*動作確認*
　`skewer -h`

<br>
<br>

**0-2-1. PRINSEQのダウンロード**
PRINSEQ(https://sourceforge.net/projects/prinseq/files/standalone/)からダウンロード

wget https://sourceforge.net/projects/prinseq/files/standalone/prinseq-lite-0.20.4.tar.gz

*解凍*
　`tar -xzvf prinseq-lite-0.20.4.tar.gz`

*prinseqディレクトリに移動*
　`cd prinseq-lite-0.20.4/`

*パーミッションの変更*
　`chmod +x *.pl`
<br>

**0-2-2. PRINSEQのPATHを通す**
`.bash_profile` に vi,vim,echo等で インストールしたPRINSEQの場所を加え、変更を反映
*（例）*
　`PATH=~/ProgramFiles/prinseq-lite-0.20.4:$PATH`
*変更したbashprofileを反映*
　`source ~/.bashrc`
*動作確認*
　`prinseq-lite.pl -h`

<br>
<br>

**0-3-1. FastQCのPATHを通す**
FastQCはサーバーにインストールされているのでPATHを通すだけでよい
`.bash_profile` に vi,vim,echo等で FastQCの場所を加え、変更を反映
　`PATH=$PATH:$HOME/.local/bin`
*変更したbashprofileを反映*
　`source ~/.bashrc`
*動作確認*
　`fastqc -h`
<br>
<br>

<br>

### 1. アダプタートリミング
skewer (参考：http://kazumaxneo.hatenablog.com/entry/2017/12/19/125341)

* single-end の場合

skewer -m tail Read.fastq.gz -o FILENAME

* paired-end の場合

skewer -m pe Read1.fastq.gz Read2.fastq.gz -o FILENAME


他にトリミングしたいadapter配列がある時は、`-x` `-y` などのオプションで除く。

FILENAME-trimmed-pair1.fastq
FILENAME-trimmed-pair2.fastq
FILENAME-trimmed.log
の3つのファイルができる。

<br>

### 2. クオリティコントロール
PRINSEQ (参考：https://bi.biopapyrus.jp/rnaseq/qc/prinseq.html)

リードの長さが 50bp未満、クオリティが 20未満の塩基をpaired-endで取り除く場合

prinseq-lite.pl -min_len 50 -trim_qual_right 20 -fastq FILENAME-trimmed-pair1.fastq -fastq2 FILENAME-trimmed-pair2.fastq -out_good FILENAME2

FILENAME2_1.fastq
FILENAME2_2.fastq
の2つのファイルに実行結果が保存される。

<br>

### 3. FastQC

FastQC（参考：https://bi.biopapyrus.jp/rnaseq/qc/fastqc.html)

結果を保存するディレクトリを先に作っておく。
`mkdir FastQC_reports`

*(基本形)*

fastqc --nogroup -o ./FastQC_reports FILENAME2_1.fastq

*(例 : 4コア使用、さらにディレクトリ内の全ての.fastqファイルやる)*

fastqc --nogroup -t 4 -o FastQC_reports *.fastq

*(例 : .fastq.gz圧縮ファイルのままでもOK)*

fastqc --nogroup -t 4 -o FastQC_reports *.fastq.gz

<br>

fastqc.html と fastqc.zip がディレクトリに保存される。
fastqc.htmlは、ブラウザで結果を確認できる。

#### 評価の例
[Good Illumina Data](http://www.bioinformatics.babraham.ac.uk/projects/fastqc/good_sequence_short_fastqc.html)
[Bad Illumina Data](http://www.bioinformatics.babraham.ac.uk/projects/fastqc/bad_sequence_fastqc.html)

<br>

### 4. Mapping
bowtie2（参考：[https://bi.biopapyrus.jp/rnaseq/mapping/bowtie2/](https://bi.biopapyrus.jp/rnaseq/mapping/bowtie2/)）

bowtie2でマッピング

bowtie2 -t -p 4 -x ~/indexのパス -1 Read1 -2 Read2 -S FILENAME.sam 2> FILENAME.summary_log.txt

samtoolsで.samを.bamに変換

samtools view -bS FILENAME.sam > FILENAME.bam

samtoolsで変換した.bamをソート

samtools sort FILENAME.bam FILENAME.sorted

FILENAME.summary_log.txtの中にマッピング率等のログが格納される。
FILENAME.sam、FILENAME.bam、FILENAME.sorted.bam ができる。

<br>
<br>
<br>
<br>
<br>
<br>
<br>
<br>
<br>
<br>
<br>
<br>
<br>
<br>

おまけ？

全部まとめシェルスクリプト(Skewer,PRINSEQ,FastQC)

(例：paired-end、標準のアダプタートリム、リードの長さが50bp未満かつクオリティが20未満の塩基を取り除く場合いてFastQCまでやる場合)

#### 1. .sh ファイルを作る。

`touch script.sh` （名前はなんでもいい）

下記をコピペ (vi, vim, echoコマンドなどで `script.sh` の中身書く)

#!/bin/bash
#$ -S /bin/bash
#$ -cwd


read1=$1
read2=$2
filename=$3

#_Skewer_adapter_trimmer_#
skewer -m pe $read1 $read2 -o $filename


#_PRINSEQ_#
prinseq-lite.pl -min_len 50 -trim_qual_right 20 -fastq $filename-trimmed-pair1.fastq -fastq2 $filename-trimmed-pair2.fastq -out_good $filename_qual


#_FascQC_#
if [ ! -d FastQC_reports ]
then
 mkdir FastQC_reports
fi 
fastqc --nogroup -t 8 -o FastQC_reports *qual_1.fastq
fastqc --nogroup -t 8 -o FastQC_reports *qual_2.fastq

パーミッションの変更
`chmod +x script.sh`

#### 2. いざq投げる！
(QNAMEと OUTNAMEは お好きなように...)

qsub -q bigmem.q -N QNAME script.sh Read1.fastq Read2.fastq OUTNAME

#### 3. あとは寝て待つ

OUTNAME.eJOB-ID、OUTNAME.oJOB-ID にエラーやらログやらが格納される