[bwa_aln_core] convert to sequence coordinate... 0.41 sec
[bwa_aln_core] refine gapped alignments... 0.07 sec
[bwa_aln_core] print alignments... 0.09 sec
[bwa_aln_core] 33353877 sequences have been processed.
[main] Version: 0.7.17-r1188
[main] CMD: bwa samse hg19.fa OUTdir15/index-1_2.valid.fastq.sai OUTdir15/index-1_2.valid.fastq
[main] Real time: 266.367 sec; CPU: 265.635 sec
[06-08 19:07:55] Running Step 3: Build Bedpe ...
buildBedpe OUTdir15/index-1_1.valid.sam OUTdir15/index-1_2.valid.sam OUTdir15/index-1.bedpe 30 1 1
Running buildbedpe...
Processed 2000000 PET
Processed 2000000 PET
Processed 2000000 PET
Processed 2000000 PET
Processed 2000000 PET
Processed 2000000 PET
Processed 2000000 PET
Processed 2000000 PET
Processed 2000000 PET
Processed 2000000 PET
Processed 2000000 PET
Processed 2000000 PET
Processed 2000000 PET
Processed 2000000 PET
Processed 2000000 PET
Processed 2000000 PET
Processed 1353876 PET
All 33353876
Reads1 Low MAPQ 22004694 65.97%
Reads1 Unmapped 8067789 24.19%
Reads2 Low MAPQ 22752786 68.22%
Reads2 Unmapped 7406492 22.21%
Output PETs 1980267 5.937%
removeDup OUTdir15/index-1.bedpe OUTdir15/index-1.rmdup.bedpe 1
Running removeDups...
Duplicates 360379
Uniques 1619888
Intra 433597
Inter 29182
OneEndMapped 1157109
buildTagAlign OUTdir15/index-1.rmdup.bedpe OUTdir15/index-1.rmdup.bedpe.tag
Running buildTagAlign...
[06-08 19:13:56] Running Step 4: MACS2 ...
macs2_wrap OUTdir15/index-1.rmdup.bedpe.tag OUTdir15/index-1 "-q 0.05"
Running Macs2...
macs2 callpeak -t OUTdir15/index-1.rmdup.bedpe.tag -f BED --keep-dup all -n OUTdir15/index-1 -q 0.05
INFO @ Tue, 08 Jun 2021 19:13:56:
# Command line: callpeak -t OUTdir15/index-1.rmdup.bedpe.tag -f BED --keep-dup all -n OUTdir15/index-1 -q 0.05
# ARGUMENTS LIST:
# name = OUTdir15/index-1
# format = BED
# ChIP-seq file = ['OUTdir15/index-1.rmdup.bedpe.tag']
# control file = None
# effective genome size = 2.70e+09
# band width = 300
# model fold = [5, 50]
# qvalue cutoff = 5.00e-02
# The maximum gap between significant sites is assigned as the read length/tag size.
# The minimum length of peaks is assigned as the predicted fragment length "d".
# Larger dataset will be scaled towards smaller dataset.
# Range for calculating regional lambda is: 10000 bps
# Broad region calling is off
# Paired-End mode is off
INFO @ Tue, 08 Jun 2021 19:13:56: #1 read tag files...
INFO @ Tue, 08 Jun 2021 19:13:56: #1 read treatment tags...
INFO @ Tue, 08 Jun 2021 19:13:57: 1000000
INFO @ Tue, 08 Jun 2021 19:13:58: 2000000
INFO @ Tue, 08 Jun 2021 19:13:58: #1 tag size is determined as 20 bps
INFO @ Tue, 08 Jun 2021 19:13:58: #1 tag size = 20.0
INFO @ Tue, 08 Jun 2021 19:13:58: #1 total tags in treatment: 2082667
INFO @ Tue, 08 Jun 2021 19:13:58: #1 finished!
INFO @ Tue, 08 Jun 2021 19:13:58: #2 Build Peak Model...
INFO @ Tue, 08 Jun 2021 19:13:58: #2 looking for paired plus/minus strand peaks...
INFO @ Tue, 08 Jun 2021 19:14:00: #2 number of paired peaks: 83978
INFO @ Tue, 08 Jun 2021 19:14:00: start model_add_line...
INFO @ Tue, 08 Jun 2021 19:14:00: start X-correlation...
INFO @ Tue, 08 Jun 2021 19:14:00: end of X-cor
INFO @ Tue, 08 Jun 2021 19:14:00: #2 finished!
INFO @ Tue, 08 Jun 2021 19:14:00: #2 predicted fragment length is 294 bps
INFO @ Tue, 08 Jun 2021 19:14:00: #2 alternative fragment length(s) may be 294 bps
INFO @ Tue, 08 Jun 2021 19:14:00: #2.2 Generate R script for model : OUTdir15/index-1_model.r
INFO @ Tue, 08 Jun 2021 19:14:00: #3 Call peaks...
INFO @ Tue, 08 Jun 2021 19:14:00: #3 Pre-compute pvalue-qvalue table...
INFO @ Tue, 08 Jun 2021 19:14:07: #3 Call peaks for each chromosome...
INFO @ Tue, 08 Jun 2021 19:14:09: #4 Write output xls file... OUTdir15/index-1_peaks.xls
INFO @ Tue, 08 Jun 2021 19:14:09: #4 Write peak in narrowPeak format file... OUTdir15/index-1_peaks.narrowPeak
INFO @ Tue, 08 Jun 2021 19:14:09: #4 Write summits bed file... OUTdir15/index-1_summits.bed
INFO @ Tue, 08 Jun 2021 19:14:09: Done!
rm OUTdir15/index-1_model.r OUTdir15/index-1_peaks.xls
[06-08 19:14:10] Running Step 5: Detect Interactions ...
extendpeak OUTdir15/index-1_peaks.narrowPeak 500 OUTdir15/index-1_peaks.slopPeak
Running extendpeak...
bedtools slop -i OUTdir15/index-1_peaks.narrowPeak -g human.hg19.genome -b 500 | bedtools merge -d 256 > OUTdir15/index-1_peaks.slopPeak
awk 'BEGIN{OFS="\t";i=1}{print $1,$2,$3,"peak_"i;i=i+1}' OUTdir15/index-1_peaks.slopPeak > tmp.bed; mv tmp.bed OUTdir15/index-1_peaks.slopPeak
peak_depth2 OUTdir15/index-1.rmdup.bedpe.tag 100 OUTdir15/index-1_peaks.slopPeak OUTdir15/index-1.peaks.depth
Running peakdepth...
psort OUTdir15/index-1.rmdup.bedpe.tag OUTdir15/index-1.rmdup.bedpe.tag.sorted > /dev/null 2>&1
bedtools coverage -sorted -b OUTdir15/index-1.rmdup.bedpe.tag.sorted -a OUTdir15/index-1_peaks.slopPeak -counts > OUTdir15/index-1.peaks.depth
build_interaction OUTdir15/index-1.rmdup.bedpe OUTdir15/index-1.peaks.depth OUTdir15/index-1.interactions OUTdir15/index-1.bedpe.stat
Running build_interaction...
pairToBed -a OUTdir15/index-1.rmdup.bedpe -b OUTdir15/index-1.peaks.depth -type both > OUTdir15/index-1.rmdup.bedpe.tmp
bedpe2Interaction OUTdir15/index-1.rmdup.bedpe.tmp OUTdir15/index-1.interactions OUTdir15/index-1.bedpe.stat
Running bedpe2Interaction...
PETs in the same peak: 24587
Intra PETs bewteen two peaks: 4065
Inter PETs bewteen two peaks: 195
Intra loops: 2246
Inter loops: 165
[06-08 19:14:28] Running Step 6: QCplots ...
Rscript /root/bin/ChIA-PET2_0.9.3/bin/QCplots.R OUTdir15 index-1
載入需要的程輯包:ggplot2
Error: package or namespace load failed for ‘ggplot2’:
package ‘scales’ was installed by an R version with different internals; it needs to be reinstalled for use with this R version
載入需要的程輯包:scales
Error: package or namespace load failed for ‘scales’:
package ‘scales’ was installed by an R version with different internals; it needs to be reinstalled for use with this R version
Error in theme_minimal() : 沒有"theme_minimal"這個函數
停止執行
(base) [[email protected] chia-pet1]# R
R version 3.6.0 (2019-04-26) -- "Planting of a Tree"
Copyright (C) 2019 The R Foundation for Statistical Computing
Platform: x86_64-redhat-linux-gnu (64-bit)
R是自由軟體,不帶任何擔保。
在某些條件下你可以将其自由散布。
用'license()'或'licence()'來看散布的詳細條件。
目标是四個結果
![](https://img.laitimes.com/img/9ZDMuAjOiMmIsIjOiQnIsIyZuBnL5MzN3EjN1UTM5AjNwEjMwIzLc52YucWbp5GZzNmLn9Gbi1yZtl2Lc9CX6MHc0RHaiojIsJye.png)
得到的前兩個想要的結果
組内的(十一列)
組間的
prefixname.interactions.MICC檔案有11+2列,最後兩列(PP( -log10(1-PostProbability) )和 FDR )是利用MICC基于貝葉斯混合模型得評估結果。