collapse taxa(example genus)

qiime taxa collapse
–i-table table.qza
–i-taxonomy taxonomy.qza
–p-level 6
–o-collapsed-table table-collapsed-L6.qza

check the number of features per taxa:

qiime feature-table summarize
–i-table table-collapsed-L6.qza
–o-visualization table-collapsed-L6-summarize.qzv
–m-sample-metadata-file sample-meta.tsv


now we do diversity core-metrics (do not miss with core-metrics-phylogenetic)

qiime diversity core-metrics \ --i-table table-collapsed-L6.qza \ --p-sampling-depth 1125 \ --m-metadata-file sample-meta.tsv \ --output-dir core-metrics_no-phylogenetics

and now we can look on PCoA:

'Bioinformatics' 카테고리의 다른 글

Convert Tsv to Biom File Format  (0) 2019.10.11
bbmap 세팅 메뉴얼  (0) 2019.01.07
NCBI blast+ local install OS Linux  (0) 2018.12.13
bedtools coverage  (0) 2018.12.04
Miso Analysis  (0) 2018.10.29

MetaGenome 분석 파이프라인에서 OTU Table을 Biom 파일 Format으로 변환하여 메모리를 줄이는 방법이 이용되고 있습니다.

 

Biom 포맷은 v1.0.0과 v2.1.0이 사용되고 있는데 이 안내는 v1.0.0에 대한 설명입니다.

 

tsv파일을 biom으로 변환하는 방법은 json형식과 hdf5형식으로 biom으로 변환이 가능합니다.

 

Usage

biom convert -i feature-table.txt -o feature-table.biom --table-type="OTU table" --to-json

biom convert -i feature-table.txt -o feature-table.biom --table-type="OTU table" --to-hdf5

 

반대의 경우. 즉 , Biom파일을 사람이 읽기 쉽게 Tsv 파일로 변환하는 방법은 다음과 같습니다.

 

Usage

biom convert -i feature-table.biom -o feature-table.txt --to-tsv

 

추가정보)

만들어진 biom 파일을 qiime2로 Import (해당 정보는 qiim2 버전에 따라 달라질 수 있음)

 

구 버전)

qiime tools import --input-path feature-table.biom --type 'FeatureTable[Frequency]' --source-format BIOMV100Format --output-path table.qza

 

확인한 바로는 --source-format BIOMV100Format은 최신 버전에서는 넣는다면 에러가 발생합니다. 따라서 다음과 같이 command를 입력하여야 합니다.

 

2018년 이후 버전)

qiime tools import --input-path feature-table.biom --type 'FeatureTable[Frequency]' --output-path table.qza

 

 

 

 

 

 

 

'Bioinformatics' 카테고리의 다른 글

PCoA- taxonomic level in Qiime2  (0) 2020.02.25
bbmap 세팅 메뉴얼  (0) 2019.01.07
NCBI blast+ local install OS Linux  (0) 2018.12.13
bedtools coverage  (0) 2018.12.04
Miso Analysis  (0) 2018.10.29

Chip seq 분석을 마친 후 narrowPeak 파일을 이용하여 Peak Profiling 을 진행하고 figure를 제작할 수 있습니다.


사용한 library는 ChIPseeker와 ggplot2 입니다.


narrowPeak 파일 format은 다음과 같습니다.


  1. chrom - Name of the chromosome (or contig, scaffold, etc.).
  2. chromStart - The starting position of the feature in the chromosome or scaffold. The first base in a chromosome is numbered 0.
  3. chromEnd - The ending position of the feature in the chromosome or scaffold. The chromEnd base is not included in the display of the feature. For example, the first 100 bases of a chromosome are defined aschromStart=0, chromEnd=100, and span the bases numbered 0-99.
  4. name - Name given to a region (preferably unique). Use '.' if no name is assigned.
  5. score - Indicates how dark the peak will be displayed in the browser (0-1000). If all scores were '0' when the data were submitted to the DCC, the DCC assigned scores 1-1000 based on signal value. Ideally the average signalValue per base spread is between 100-1000.
  6. strand - +/- to denote strand or orientation (whenever applicable). Use '.' if no orientation is assigned.
  7. signalValue - Measurement of overall (usually, average) enrichment for the region.
  8. pValue - Measurement of statistical significance (-log10). Use -1 if no pValue is assigned.
  9. qValue - Measurement of statistical significance using false discovery rate (-log10). Use -1 if no qValue is assigned.
  10. peak - Point-source called for this peak; 0-based offset from chromStart. Use -1 if no point-source called.

1~6 column은 Bedfile format과 동일 합니다. 

track type=narrowPeak visibility=3 db=hg19 name="nPk" description="ENCODE narrowPeak Example"
browser position chr1:9356000-9365000
chr1    9356548 9356648 .       0       .       182     5.0945  -1  50
chr1    9358722 9358822 .       0       .       91      4.6052  -1  40
chr1    9361082 9361182 .       0       .       182     9.2103  -1  75



아래 코드들은 R 에서 작성되었습니다.


Usage1


library(ChIPseeker)

library(ggplot2)

peak <- readPeakFile("narrowPeak파일")
peak
## GRanges object with 1331 ranges and 2 metadata columns:
##          seqnames              ranges strand |             V4        V5
##             <Rle>           <IRanges>  <Rle> |       <factor> <numeric>
##      [1]     chr1       815093-817883      * |    MACS_peak_1    295.76
##      [2]     chr1     1243288-1244338      * |    MACS_peak_2     63.19
##      [3]     chr1     2979977-2981228      * |    MACS_peak_3    100.16
##      [4]     chr1     3566182-3567876      * |    MACS_peak_4    558.89
##      [5]     chr1     3816546-3818111      * |    MACS_peak_5     57.57
##      ...      ...                 ...    ... .            ...       ...
##   [1327]     chrX 135244783-135245821      * | MACS_peak_1327     55.54
##   [1328]     chrX 139171964-139173506      * | MACS_peak_1328    270.19
##   [1329]     chrX 139583954-139586126      * | MACS_peak_1329    918.73
##   [1330]     chrX 139592002-139593238      * | MACS_peak_1330    210.88
##   [1331]     chrY   13845134-13845777      * | MACS_peak_1331     58.39
##   -------
##   seqinfo: 24 sequences from an unspecified genome; no seqlengths


출처 : https://genome.ucsc.edu/FAQ/FAQformat.html#format12 





Usage2


peak <- readPeakFile("narrowPeak파일")
covplot("narrowPeak파일")



♣각각의 chr에서의 peak을 profile 할 수 있습니다.




Usage3



peak = GenomicRanges::GRangesList(샘플명=readPeakFile("narrowPeak파일"), 샘플명=readPeakFile("narrowPeak파일")
covplot(peak, weightCol="V5") 

 

♣Control과 Test의 peak의 차이를 한번에 확인할 수 있습니다.



Usage4


peak = GenomicRanges::GRangesList(샘플명=readPeakFile("narrowPeak파일"), 샘플명=readPeakFile("narrowPeak파일")
covplot(peak, weightCol="V5", chrs = c("chr선택","chr선택"), xlim=c(x축 start position, X축 end position)) 



♣ 원하는 각 chromosome의 size를 제한 하여 peak을 볼 수 있습니다.



Additional Options


xlab  -> x축 title

ylab -> y축 title

title -> figure title

----------------------------------


### + facet_grid(chr ~.id)  -> 범례를 설정합니다. 또한 y축의 scale이 setting 됩니다.



'Bioinformatics > Chip-seq' 카테고리의 다른 글

Heatmap을 이용한 Chip-seq data visualization  (0) 2018.10.11
Chip-seq 분석이란?  (0) 2018.10.11

+ Recent posts