(Cut and Mark / Cut and Run) 잘라내기 및 표시 튜토리얼

by

기사 https://yezhengstat.github.io/CUTTag_tutorial/#Cite_this_tutorial 인용하다.

CUT & TAG 데이터 처리 및 분석 자습서

1. 소개

진핵 핵 DNA에서 발생하는 모든 동적 프로세스는 염색질 지형/토폴로지를 따릅니다.

염색질 환경은 뉴클레오솜, 그 변형, 전사 인자 및 염색질 관련 복합체로 구성됩니다.

이러한 독특한 염색질 기능을 사용하여 다양한 세포 유형과 발달 단계 사이의 활성화/침묵 부위 또는 염색질 도메인을 표시할 수 있습니다.


염색질 도메인

절단 및 표시 방법


  1. 세포 또는 핵을 투과화(항체 통과 허용)
  2. 항체와 함께 배양 → 특정 표적에 결합.
  3. pA-Tn5+ 어댑터는 항체 결합 부위에 연결(도킹)됩니다.

  4. 마그네슘 철을 사용하여 트랜스포솜이 활성화되고 트랜스포솜의 어댑터가 주변 DNA에 통합됩니다.

  5. DNA를 추출하고 증폭(PCR)했습니다.

장점: 항체 테더링된 Tn5 기반 방법은 pA-Tn5 테더링 후 샘플의 엄격한 세척 및 어댑터 통합의 고효율로 인해 높은 감도를 달성합니다.

(Chipseq의 신호 대 잡음비를 개선할 수 있음)

→ 최대 90개 샘플 가능

지도 시간

인간 림프종 K562 세포주에서 히스톤 변형 프로파일링.

** 염색질 단백질(TF, RNA polymerase II, epitope tag protein)과 관련된 모든 분석이 가능합니다.


# 필요한 packages 
library(dplyr)
library(stringr)
library(ggplot2)
library(viridis) # https://m.blog.naver.com/regenesis90/222201170325
library(GenomicRanges)
library(chromVAR) ## For FRiP analysis and differential analysis
library(DESeq2) ## For differential analysis section
library(ggpubr) ## For customizing figures
library(corrplot) ## For correlation plot

2. 데이터 전처리

QC

  • 빠른 QC
  • “염기 서열 내용당” – CUT&TAG 실험을 수행할 때 처음 10bp에 들쭉날쭉한 불일치가 포함될 수 있습니다.

    이는 정렬 또는 피크 호출을 수행할 때 문제가 되지 않습니다.

    (튜토리얼 Bowtie2 매개변수를 사용하면 정확한 매핑 정보를 제공하므로 트리밍을 권장하지 않습니다.

    )

(필요한 경우 기술 복제/레인 병합)

3. 정렬

단일 HiSeq 2500 플로우 셀에서 최대 90개의 풀링된 샘플에 대한 단일 인덱스 25×25 PE Illumina 시퀀싱(끝당 25bp?)을 수행할 수 있습니다.

각 샘플에는 고유한 PCR 프라이머 바코드가 있습니다.

각 라이브러리가 약 500만 페어 엔드 읽기를 갖도록 조정하여 고품질 분석을 달성할 수 있습니다.

나비 넥타이 2 정렬


https://bowtie-bio.sourceforge.net/bowtie2/manual.shtml#the-bowtie2-aligner

##== linux command ==##
cores=8
ref="/path/to/bowtie2Index/hg38"

## Build the bowtie2 reference genome index if needed:
## bowtie2-build path/to/hg38/fasta/hg38.fa /path/to/bowtie2Index/hg38©

bowtie2 --end-to-end --very-sensitive --no-mixed --no-discordant --phred33 -I 10 -X 700 -p ${cores} -x ${ref} -1 ${projPath}/fastq/${histName}_R1.fastq.gz -2 ${projPath}/fastq/${histName}_R2.fastq.gz -S ${projPath}/alignment/sam/${histName}_bowtie2.sam &> ${projPath}/alignment/sam/bowtie2_summary/${histName}_bowtie2.txt

-end-to-end –매우 민감한 –no-mixed –no-discordant –phred33 -I 10 -X 700:10-700bp 길이 매핑을 삽입합니다.

25×25 PE 시퀀싱을 하는 경우 트리밍할 필요가 없지만 더 긴 시퀀싱을 하는 경우 Cutadapt 또는 -local –매우 민감함 –no-mixed –no-discordant –phred33 -I 10 -X 700 사용을 위해서는 읽기의 3′ 끝에 남아 있는 어댑터를 제거해야 합니다.

(-local을 사용하는 경우 Bowtie2는 정렬 점수를 최대화하기 위해 끝을 자를 수 있습니다.

)

스파이크인 교정을 위한 스파이크인 게놈과의 정렬

반응 동안 대장균 DNA는 대장균 DNA를 포함하는 박테리아에 의해 생성된 pA-Tn5 단백질과 함께 운반되며 무작위로 표지됩니다.

…하지만 효과가 좋지 않아 아직은 소용이 없다는 분들이 많습니다.

비교 요약 bowite2 매뉴얼

2984630 reads; of these: # sequencing depth (total number of paired reads)
  2984630 (100.00%) were paired; of these:
    125110 (4.19%) aligned concordantly 0 times
    2360430 (79.09%) aligned concordantly exactly 1 time
    499090 (16.72%) aligned concordantly >1 times # 2360430 + 499090 = successfully mapped reads
95.81% overall alignment rate
  • 80%, 100만 개의 매핑된 단편은 인간 게놈의 히스톤 변형에 대한 신뢰할 수 있는 개요를 제공합니다.


중복 제거

pA-Tn5가 링커를 DNA에 삽입할 때 정확한 삽입 위치는 주변 DNA의 접근성에 영향을 받습니다.

따라서 정확히 동일한 위치에 삽입되어 동일한 조각이 생성되며 이를 “중복”이라고 합니다.

(이는 PCR에 의해 생성된 중복과는 다릅니다.

) 그러나 고품질 CUT&RUN 데이터셋의 중복률은 매우 낮고 추가된 ‘중복’ 조각은 실제 조각이므로 잘 제거되지 않습니다.

그러나 데이터가 매우 작거나 PCR에서 중복이 의심되는 경우 제거해야 합니다.

# example markedDup.bam
$samtools view ./A1_H3K27_IRI_markedDup.bam | head -1
K00124:539:H7NVGBBXX:7:1110:25083:12314 99      scaffold_1      13      11      1S94M1S =       13      96      AACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCGTTCCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCT        AAFFFJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJ MD:Z:17T0T75    PG:Z:MarkDuplicates     XG:i:0  NM:i:2  XM:i:2  XN:i:0  XO:i:0  AS:i:172        XS:i:174        YS:i:173        YT:Z:CP

**고유 조각 수 = MappedFragNum_hg38 * (1-DuplicationRate/100)


맵 조각 크기 분포 평가

CUt & RUN 라벨링 실험의 경우 조각 크기는 뉴클레오솜 크기의 조각(~180bp 또는 x2, x3)이거나 더 짧아야 합니다(인자 결합 부위).

**정렬된 sam 파일의 9열에서 조각 길이를 확인할 수 있습니다.

## Extract the 9th column from the alignment sam file which is the fragment length
samtools view -F 0x04 $projPath/alignment/sam/${histName}_bowtie2.sam | awk -F'\\t' 'function abs(x){return ((x < 0.0) ? -x : x)} {print abs($9)}' | sort | uniq -c | awk -v OFS="\\t" '{print $2, $1/2}' >$projPath/alignment/sam/fragmentLen/${histName}_fragmentLen.txt

4. 비교 결과 필터링 및 파일 형식 변환

도면 품질

samtools view -q minQualityScore를 사용하여 minQualityScore보다 낮은 경우 제거할 수 있습니다.

이 시점에서 QualityScore는 MAPQ(x) ****= -10 * log_10(p)입니다.

(p: x가 잘못 매핑될 확률)

파일 형식 변환

## Filter and keep the mapped read pairs
**samtools view -bS** -F 0x04 $projPath/alignment/sam/${histName}_bowtie2.sam >$projPath/alignment/bam/${histName}_bowtie2.mapped.bam
	# -F 0x04: exclude 0x04 flag = reads that are unmapped for the output 

## Convert into bed file format
**bedtools bamtobed** -i $projPath/alignment/bam/${histName}_bowtie2.mapped.bam -bedpe > $projPath/alignment/bed/${histName}_bowtie2.bed
	# example bed line 
	# scaffold_16     71719668        71719805        scaffold_16     71719676        71719825        K00124:539:H7NVGBBXX:7:1101:1225:64966  44      +       -

## Keep the read pairs on the same chromosome and fragment length < 1000bp.
awk '$1==$4 && $6-$2 < 1000 {print $0}' $projPath/alignment/bed/${histName}_bowtie2.bed > $projPath/alignment/bed/${histName}_bowtie2.clean.bed
	# awk: for text processing and pattern matching
		# $1==$4: 1st field == 4th field 인지 확인 (same chromosome)
		# $6 position = second pair's start position 
		# &&(and)
	# If both conditions are true, the awk script outputs the entire record ($0), which consists of all the fields in the current line of the input file.

## Only extract the fragment related columns
cut -f 1,2,6 $projPath/alignment/bed/${histName}_bowtie2.clean.bed | sort -k1,1 -k2,2n -k3,3n  > $projPath/alignment/bed/${histName}_bowtie2.fragments.bed
	# first field (-k1,1) in ascending order (-k2,2n) and then based on the second field (-k2,2n) and third field (-k3,3n) in numeric order

액세스 복제 재현성

복제 사이의 재현성을 조사하기 위해 게놈을 500bp 빈으로 나누고 각 복제에 대해 계산되고 비교되는 log2 변환된 읽기 수의 Pearson 상관관계를 나타낼 수 있습니다.

## We use the midpoint of each fragment to infer which 500bp bins this fragment belongs to.
binLen=500
awk -v w=$binLen '{print $1, int(($2 + $3)/(2*w))*w + w/2}' $projPath/alignment/bed/${histName}_bowtie2.fragments.bed | sort -k1,1V -k2,2n | uniq -c | awk -v OFS="\\t" '{print $2, $3, $1}' |  sort -k1,1V -k2,2n  >$projPath/alignment/bed/${histName}_bowtie2.fragmentsCount.bin$binLen.bed


6. 피크 콜링

SEAC

Cut & Run(SEART) 도구를 사용한 희소 농축 분석을 사용하여 염색질 프로파일링 데이터에서 호출 피크 및 농축 영역을 식별할 수 있습니다.

SEACR은 두 가지 유형의 피크 통화에 사용할 수 있습니다.

(1) 좁은 피크 – 인자 결합 부위의 경우 (2) 넓은 피크 – 히스톤 수정 호출의 경우

-입력하다: 침대 차트 파일(페어드 엔드 시퀀싱에서). 면역글로불린 침대 지도 파일

– Do: IgG 컨트롤 데이터와 비교하여 피크 높은 신호 강도 피크.

SEAR 매개변수:

– 대장균 읽기 횟수(스파이크 인)를 사용하여 조각을 정규화하는 경우 “not”

– 정규화된 “규범”

seacr="SEACR/SEACR_1.3.sh"
histControl=$2
mkdir -p $projPath/peakCalling/SEACR

bash $seacr $projPath/alignment/bedgraph/${histName}_bowtie2.fragments.normalized.bedgraph \\
     $projPath/alignment/bedgraph/${histControl}_bowtie2.fragments.normalized.bedgraph \\
     non stringent $projPath/peakCalling/SEACR/${histName}_seacr_control.peaks

bash $seacr $projPath/alignment/bedgraph/${histName}_bowtie2.fragments.normalized.bedgraph 0.01 non stringent $projPath/peakCalling/SEACR/${histName}_seacr_top0.01.peaks

호출된 봉우리의 수

복제당 피크 수를 구성하고 조사하면 데이터를 이해하는 데 도움이 됩니다.

→ 생물학적 복제(rep1, rep2 중첩 피크)의 재현성으로 데이터의 신뢰도를 확인할 수 있습니다.

피크 영역(FRiP)의 조각 비율

FRiP는 신호 대 잡음비의 척도입니다.

20회 이상 통과해야 합니다.

피크 수, 피크 폭, 피크 재현성 및 FRiP의 시각화