Chromatin Immunoprecipitation Sequencing Analysis Report

Novogene Co.,Ltd.

I.Background

Chromatin Immunoprecipitation Sequencing (ChIP-Seq) provides genome-wide profiling of DNA targets for histone modifications, transcription factors, and other DNA-associated proteins. It combines the selectivity of chromatin immuno-precipitation (ChIP) for recovering specific protein-DNA complexes, with the power of next-generation sequencing (NGS) for high-throughput sequencing of the recovered DNA. Additionally, because the protein-DNA complexes are recovered from living cells, binding sites can be compared in different cell types and tissues, or under different conditions. Applications range from transcriptional regulation to developmental pathways to disease mechanisms and beyond.

Novogene experienced scientists offers comprehensive analysis solutions to meet your project needs by ensuring data and analysis quality that is publication ready.

Novogene Co.,Ltd.

II. Experimental Procedure

From DNA samples to final data report, each step of sample QC, library preparation and sequencing may have an impact on the quality and quantity of data output, and the quality of data will directly affect the results of subsequent information analysis. Therefore, obtaining high-quality data is the premise to ensure that bioinformatics analysis is correct, comprehensive and credible. In order to ensure the accuracy and reliability of sequencing data from the beginning, Novogene strictly checks every experimental step and fundamentally ensures the output of high-quality data. The workflow is as follows:

Project Workflow

1　DNA sample test

The DNA samples are immunoprecipitated and purified on the end of client. Novogene detection of these IPed DNA samples mainly includes two steps:

(1) Qubit accurately quantifies DNA concentration, based on which the total amount was calculated.

(2) Agarose gel electrophoresis analyses of DNA fragments size distribution and detects any contamination.

2　Library Preparation for Sequencing

After the DNA samples quality control:

(1) The DNA fragments are repaired, dA-tailed.

(2) The DNA fragments with A tail are ligated to sequencing adaptors.

(3) The final DNA library is obtained by size selection and PCR amplification.

The workflow of library construction is as follows:

Library Construction Workflow

3　Quality control of Library

After the construction of the library, the initial quantification was done with Qubit 2.0, and the library is diluted to 1ng/l. Then the insertion size of the library is detected with NGS3K. If meeting the expectation, the accurate concentration of the library is quantified by Q-PCR (library effective concentration > 2nM) to ensure the accurate molar amount that will pooled for sequencing.

4　Sequencing

After library quality control, sequencing is performed for different libraries according to the concentration and the demand of data amount on Illumina NovaSeq platform. The basic principle of sequencing is sequencing by synthesis. Four kinds of fluorescent labeled dNTP, DNA polymerase and adapter primers are added to the flow cell for amplification. When each sequence cluster extends the complementary chain, each fluorescent labeled dNTP releases the corresponding fluorescence. The sequencer captures the fluorescent signal and converts the light signal into the sequencing peak through computer software to obtain the sequence information of the fragment to be detected. The sequencing process is shown in the following figure.

The sequencing principle and process

Novogene Co.,Ltd.

1.3　FastQC of Clean_reads

Quality control of the clean reads after trimming is in the following:

si_Ctrl_in si_Ctrl si_SIRT3 si_SIRT3_in

Figure 1-2 FastQC of Clean_reads

Top-left: horizontal axis represents position in read (bp), vertical axis represents quality scores, different colors represent different quality range.

Top-middle: horizontal axis represents position in read (bp), vertical axis represents base percentage.

Top-right: the distribution of read length.

Bottom-left: mean GC content (%). Red line is the actual GC content distribution, blue line is the theoretical distribution.

Bottom-middle: N-ratio across all bases.

Bottom-right: x-axis is the sequence duplication level. y-axis is the percent of the deduplicated reads.

Novogene Co.,Ltd.

2　Mapping

2.1　Summary of mapping

The common tools for mapping are Bowtie, BWA, MAQ, TOPhat, etc. We choose proper softwares and parameters according to different genome characters to do the genome mapping analysis for the filtered reads. Considering the small fragment size of ChIP-Seq, and the percentage of the unique sequence in the total sequence is the most important information, thus, we can map the reads to the reference genome using BWA much more accurately (Li, H. and R. Durbin, 2009). The summary of mapping is in the following table:

Table 2-1 Summary of mapping

Unique_mapped: reads with MAPQ (Li and Ruan et al., 2008) no lower than 13; can be interpreted as the chance of non-accurate mapping (same score for the random mapping) is 0.05.

Duplicates: the reads mapped to the exact same position of the genome;

Mapped is relatively to Clean_reads;

Unique_mapped is relatively to Mapped;

Dup_Unique_mapped is relatively to Unique_mapped.

2.2　MAPQ

The most important thing during ChIP-Seq analysis is the percentage of the unique sequence in the total sequence number. Duplicates are labeled using SAMBLASTER (Faust and Hall, 2014) and mapping quality value is calculated by MAPQ. Proper quality value is chosen as the only threshold for mapping. Here we choose 13 as the threshold, which refer to that the mapping chance of the accordingly non-unique region is only 0.05. Only one read pair is kept for the duplicates in the followed peak calling.

si_Ctrl_in si_Ctrl si_SIRT3 si_SIRT3_in

Figure 2-1 Distribution of MAPQ

Horizontal axis is MAPQ, vertical axis is the reads count.

Novogene Co.,Ltd.

2.3　Genome-wide distribution of the mapped reads

Summary of the density of total mapped reads in different chromosomes (plus and minus) is shown in the followed figure. With 5k slide window size, the medium of the number of reads mapped to each base within window is calculated, and converted to log₂. The longer the whole chromosome is, the more reads are mapped (Marquez et al. 2012). From the correlation of the number of the mapped reads and the length of the chromosome in the figure, we can see the correlation of the total number of the reads and the length of the chromosome much easier.

si_Ctrl_in si_Ctrl si_SIRT3 si_SIRT3_in

Figure 2-2 Genomewide distribution of the mapped reads

Horizontal axis represents the postition of the chromosome, vertical axis represents the number of the reads mapped to 1000nt window size. Here is the unique mapping and deduplication results.

2.4　Distribution of the reads mapped to the gene

Since the binding sites of transcription factor and histone protein are important for gene regulation, analysis of relative mapping position distribution can help us predict the protein function. Each gene and its 2kb upstream and 2kb downstream are divided into 100 equal parts. Reads density are calculated in each part by the mapped reads and the ratio of the reads to total reads.

si_Ctrl_in si_Ctrl si_SIRT3 si_SIRT3_in

Figure 2-3 Distribution of the reads mapped to gene

Horizontal axis: relative position of the gene. Vertical axis：reads density.

Novogene Co.,Ltd.

2.5　Sample correlation detection

Biological replicate is necessary for every experiment, same for high-throughput technology (Hansen et al.). Biological replicate mainly has two applications: One is to prove that the biological experiment can be replicated and there is no large variance. The other one is to make sure that following differential gene analysis can get reliable results. The correlation among samples is an important index to see whether the experiment design is reliable and whether the sampling is right. The correlation coefficient is much closer to 1, the similarity of the reads distribution pattern among samples is much higher.

Figure 2-4 Spearman test among samples

Heat-map of sample correlation test, correlation coefficient among samples (Spearman correlation coefficient). The darkness of the color represents how large is the correlation coefficient.

2.6　Visualization of pileup signal

We provide the visualization results of genome wide reads mapping in bam format. IGV (Integrative Genomics Viewer) browser is recommended to view the bam file. IGV browser has the following characters: (1) can reveal single or multiple mapping positions in the genome in different scales, including the distribution of the reads in different chromosomes, and the distribution of annotated exons, introns, splicing junctions and inter-gene region; (2) can reveal reads abundance in different region under different scales; (3) can reveal the annotation information of the gene and alternative splicing isoforms; (4) can reveal other annotation information; (5) can download annotation information from remote and local server. Please check the IGV manual for detail instruction (IGVQuickStart.pdf).

Figure 2-5 Visualization of pileup signal by IGV(demo)

Novogene Co.,Ltd.

3　Fragment size prediction

3.1　Summary of fragment size

For a specific binding site, there is a significant reads enrichment in the binding site. For single-end sequencing, we use MACS2 software to predict the frag_sizes of IP experiment. The whole genome is scaned using certain window size and the enrichment level of the reads in each window is calculated. Then proper windows (eg.1000) is extracted as the samples to build the enrichment model to predict the length of frag_sizes. For double-end sequencing, we use RSeQC software to predict the frag_sizes for the mapping results. The predicted frag_sizes are used for later peak calling.

Table 3-1 Summary of Fragment size

Sample: sample name

frag_sizes_length: frag_sizes length

infor: default parameter

3.2　K distribution of fragment size

Length distribution is analyzed by kernel density method and shown below:

si_Ctrl_in si_Ctrl si_SIRT3 si_SIRT3_in

Figure 3-1　frag_size distribution Figure

Single end:

left: horizontal axis is the distance to the center of peak model. vertical axis is the percentage of the reads chosen for modeling.

Red color represents plus strand, blue color represents minus strand.

Right: horizontal axis is the distance to the middle of the peak model. vertical axis is the correlation between plus strand and minus strand. The distance from the red dash line to the middle is the predicted frag_sizes.

Double end:

horizontal axis is the length of the predicted length of frag_sizes, vertical axis is the value of kernel density.

Novogene Co.,Ltd.

4　Strand cross correlation

4.1　Summary of strand cross correlation

As we know that the sequenced reads will be approximately evenly distributed on the plus and minus strands. By calculating the correlation of plus and minus strand (SCC), we can get the distance between two strands. Through the SCC analysis of IP and input data, we can not only obtain the correlation coefficient between plus and minus strands, but also test the effect of IP experiment.

Table 4-1 Summary of strand cross correlation

Sample: sample name

Median_read_length: mean value of reads length

Predicted_fragment_length: predicted length of fragment sizes

CC_min: the lowest SCC

CC_read_length: the SCC of the longest reads

CC_fragment_lengt: fragment sizes relative SCC

NSC: normalized strand coefficient

RSC: relative strand correlation

4.2　Plots of strand cross correlation

The normalized SCC curve for all the samples is shown below:

Figure 4-1 Plots of strand cross correlation

Horizontal axis is the length of strand shift. Vertical axis is the value of normalized strand cross correlation. Different color represents different samples.

Novogene Co.,Ltd.

4.3 SCC distribution between experimental groups

The ratio of the CC (Cross Correlation) value of the same experimental group (including IP and Input) to the lowest CC value of the whole SCC curve (NSC) should be no less than 1.05. Besides, there is another shadow peak in the reads which is resulted from the sequencing read length. RSC value refers to the ratio of two differences: the difference between the CC value of fragment length and the lowest CC value, and the difference between CC value in the shadow peak and the lowest CC value.

CTRL KD si_Ctrl_in si_Ctrl si_SIRT3 si_SIRT3_in

Figure 4-2　SCC curve

Horizontal axis is the lag between plus and minus tags during Pearson Correlation Coefficient calculation. Vertical axis is the Pearson Correlation Coefficient.

Novogene Co.,Ltd.

5　Peak Calling

The annotation of transcription factor binding sites, histone binding sites is the important information for understanding the regulation mechanism and function. By mapping to the reference genome we can obtain the information of protein-DNA binding sites directly. By making use of MACS2 software (Yong Zhang, Tao Liu et al., 2008) (threshold q value = 0.05) to finish the peak calling, we can calculate the number of peaks, the peak width and its distribution, and find the peak related genes.

5.1　Summary of peak calling

Table 5-1 Summary of peak calling

Experiment: experimental group name(one ChIP Experiment includes an IP and a control, eg. Input or Mock or no control);

EIP:experiment name after chip handling;

EInput/Mock: control group;

EFragment_length: frag sizes length

ECount_of peak: the number of peak (narrow). If would like to test broad peak, need annotate in the information collection form;

EFRiP: the ratio of the numer of the reads in the peak to the total reads, which can test the effect of IP experiment;

ECount_of_summits: number of summits. Some peaks can have multiple summits due to close position.

Novogene Co.,Ltd.

5.2　Genome wide distribution of peaks

The summary of genome wide distribution of peaks is shown in the figure below. From the number of the peaks mapping to the chromosome and its distribution can reflect the distribution of the protein binding sites. When the number of the chromosomes with peaks is larger than 30, only 15 chromosomes are shown.

CTRL KD

Figure 5-1　Genome wide distribution of peaks

Horizontal axis is the coordinate of the peaks in the chromosome. Vertical is the chromosomes. Every blue bar represents a peak.

5.3　Distribution of peak width

The peak width represents the length of the DNA that is bound by the studied protein. The number of peaks of different peak width is calculated, as shown below.

CTRL KD

Figure 5-2　Distribution of peak width

Horizontal axis is the width of the peak (nt), vertical axis is the number of the corresponding peaks.

Novogene Co.,Ltd.

5.4　Distribution of fold enrichment

The fold enrichment value is the digital display of the peak signal during peak calling, which can be called signal value. The larger the value, the more reads are enriched to this peak. The peak number of different fold changes is calculated. The distribution of fold enrichment is shown below:

CTRL KD

Figure 5-3　Distribution of fold enrichment

Horizontal axis represents the enrichment fold change of the peak. Vertical axis represents the number of the peaks.

5.5　Distribution of q values

The significance of the peak is the indication of the credibility of the peak, so the q value for eah peak is calculated. The peak number of differnet q values is calculated and the distribution is shown below:

CTRL KD

Figure 5-4　Distribution of q values

Horizontal axis represents -log₁₀ q value of the peak; Vertical axis represents the number of the peaks.

Novogene Co.,Ltd.

5.6　Count of summits in peaks

The number of summits in each peak is analysized, which infers the peak type in the IP experiment.

CTRL KD

Figure 5-5　Count of summits in peaks

Horizontal axis represents the number of summit in each peak. Vertical axis is the corresponding number of the peak.

5.7　Summits distribution

Each peak is divided to 100bp windows and the summits in each window of all the peaks are counted as shown below:

CTRL KD

Figure 5-6　Percentile position of summits in peaks

Horizontal axis is the position of sumits, vertical axis is the count of the summits.

Novogene Co.,Ltd.

6　Motif analysis

The binding of protein, such as transcription factor, histone etc. to DNA is not random, but sequence preference. MEME(Timothy L. Bailey and Charles Elkan, 1994)and Dreme (Timothy L. Bailey, 2011) softwares is used to detect significant motif sequence in the peak. Tomtom (Shobhit Gupta, JA Stamatoyannopolous, 2007) software is used to annotate the motif by mapping it to the annotated Motif database. We use sequence logo to show the base bias in different position in the binding sites in long Motif (8~30) (Fig 6-1) and short Motif (~8) (Fig 6-2). (Note: because of binding site specification, one motif sequence can be only shown in one section (<= 8 or >= 9), and parts of the figure below may have empty results).

6.1　Motif searching

CTRL KD

Figure 6-1　long conservative sequence

Logos are listed in order. The motif in the right is the reverse complementary sequence of the left one.

CTRL KD

Figure 6-2　short conservative sequence

Logos are listed in order. The motif in the right is the reverse complementary sequence of the left one.

Novogene Co.,Ltd.

6.2　MEME motif annotation

Tomtom is used to annotate the motif sequences detected by MEME. (Note: since the conservative sequence in the binding site is short, motif can be only in one region (<= 8), resulting in no result for meme test).

Table 6-1 MEME detection Motif calculation

#Query ID: detected motif;

Target ID: known motif ID in the database;

Optimal offset: the number of lag bases;

p-value:probability of MCMC;

E-value:false positive probability;

q-value: FDR value;

Overlap: overlapping base pair between two sequences;

Query consensus: detected motif sequence;

Target consensus: motif sequence in the target database;

Orientation: plus or minus strand for the target sequence.

The sequence logo below shows the comparison results between MEME detected motif and known motif. The result is in the following:

CTRL KD

Figure 6-3　MEME motif annotation

Novogene Co.,Ltd.

6.3　DREME motif annotation

DREME is a motif discovery algorithm specifically designed to find the short, core DNA-binding motifs of eukaryotic TFs. (Note: Since the conservative sequence in the binding sites is long, the motif could only show up in one region (<= 8), which could show empty result with Dreme detection)

Table 6-2 Summary of motif detected by Dreme

#Query ID: detected motif ID;

Target ID: known motif ID in the database;

Optimal offset: lag base number;

p-value: probability of MCMC;

E-value: false positive probability;

q-value: FDR value;

Overlap: overlapping base between two sequences;

Query consensus: detected motif base sequence;

Target consensus: motif base sequence in the target database;

Orientation: plus or minus tag for the targeted sequence.

The sequence logo below shows the comparison result between detected motif with Dreme and known motif:

CTRL KD

Figure 6-4　Dreme motif annotation

Novogene Co.,Ltd.

7　Peak annotation

7.1　Peak-TSS distance

Peak-TSS distance distribution can be used to predict the protein binding sites. On the one hand, the effect of IP assay can be detected according to the protein binding sites; on the other hand, the regulation mode or function of the protein can be predicted through the characteristics of protein binding. TSS (transcription start site) of every peak related gene is detected by PeakAnnotator (Salmon-Divon and Dvinge et al., 2010). Peak-TSS distance distribution is analyzed according to the distance between the peak and TSS.

CTRL CTRL_enlarged KD KD_enlarged

Figure 7-1　Distribution of Peak-TSS Distance

Horizontal axis represents the distance from peak to TSS. Vertical axis represents the number of the peak.

7.2　Distribution of peaks in functional regions

The figure below shows the distribution of peaks in different functional regions.

CTRL KD

Figure 7-2 Peak distribution in functional gene region

Distribution of peak in different functional area. Horizontal axis represents different functional area, vertical axis represents the ratio of the peak in the functional region to the total peaks. The number on the top of the functional region represents peak number.

U2000 means 2000bp in the upstream region, D2000 means 2000bp in the downstream region;

CDSu2K and CDSd2K means upstream and downstream 2kb of CDS;

TSS100, Stop100, start100 means 100bp centered with TSS,TTS, Start-codon and stop-codon.

Novogene Co.,Ltd.

7.3　GO enrichment analysis

Gene Ontology (GO, http://www.geneontology.org/) is a major bioinformatics initiative to unify the representation of gene and gene product attributes across all species. More specifically, the project aims to: 1) maintain and develop its controlled vocabulary of gene and gene product attributes; 2) annotate genes and gene products, and assimilate and disseminate annotation data; and 3 provide tools for easy access to all aspects of the data provided by the project, and to enable functional interpretation of experimental data using the GO, for example via enrichment analysis. GO analysis covers three aspects: Molecular Function, Biological Process and Cellular Component.

Results of GO enrichment analysis of the peak related gene are showed below:

Table 7-1 peak related gene GO enrichment

GO_accession: the unique GO ID;

Description: Function description;

Term_type: Including cellular_component, biological_process and molecular_function;

Over_represented_pValue: significance of enrichment analysis;

Corrected_pValue: P-Value after correction, normally, with padj< 0.05 are enriched;

Gene_with_peak_item: the number of peak related genes with this GO term;

Gene_with_peak_list: the number of peak related gene;

Bg_item: the number of background gene with this GO term;

Bg_list: the number of background genes;

Genes: Annotated peak related gene ID (not show in this table);

Choose the first five results for each experiment.

Novogene Co.,Ltd.

Below is the GO enrichment bar of peak overlapping genes. It is directly reflect the distribution of the number of the peak overlapping genes that are enriched in the biological process, cellular component and molecular function.

CTRL CTRL_mf CTRL_cc CTRL_bp KD KD_mf KD_cc KD_bp

Figure 7-3　GO enrichment

Bar: horizontal axis is the number of the peak overlapping gene. Different color is used to differentiate biological process, cellular component and molecular function.

DAG figure: every note represents a GO project. Rectangular represents top 10 enrichment GO. Darkness represents the enrichment level. The darker of the color, the higher of the enrichment.

The name of the project and the significance (padj) are shown in every note.

Novogene Co.,Ltd.

7.4　KEGG enrichment analysis

Different genes coordinate to each other to realize the function in organism. Kyoto Encyclopedia of Genes and Genomes (KEGG, http://www.kegg.jp/) is a collection of manually curated databases, containing resources onf genomic, biological-pathway and disease information (Kanehisa,2008). Pathway enrichment analysis identifies significantly enriched metabolic pathways or signal transduction pathways associated with peak overlapping genes, comparing the whole genome background.

Table 7-2 KEGG enrichment list

(1) Term: Description information of KEGG pathway.

(2) Database: KEGG PATHWAY.

(3) ID: The only pathway identifier in the KEGG database.

(4) Input number: the number of peak overlapping gene in the corresponding pathway.

(5) Background number: The number of the genes in the corresponding pathway.

(6) P-value: statistical significant level of enrichment analysis.

(7) Corrected p-value: statistical significant level after correction. Normally, if the corrected P-value < 0.05 this pathway is enriched.

(8) Input: The peak overlapping genes in the corresponding pathway. Because of too much genes, here we ignore it.

(9) Hyperlink: see result file.

Peak overlapping gene KEGG enrichment scatter plot (Figure 3.24) is the visualization of KEGG enrichment analysis. In this figure, the extent of KEGG enrichment is measured by Rich factor, qvalue and the number of the genes that are enriched in this pathway. The top 20 significantly enriched terms in the KEGG enrichment analysis are displayed below. If the enriched pathway is less than 20, all the terms will be displayed.

CTRL KD

Figure 7-4　KEGG enrichment scatter

Vertical axis represents pathway name, horizontal axis represents Rich factor, the size of the dot represents the numbers of the overlapping peak genes in the pathway, and the color of the dots corresponds to different qvalue range. The size of the dots represents the number of the genes whose peaks are overlapping in this pathway.

Novogene Co.,Ltd.

8　Differential analysis

Only when the number of groups is not less than two can the difference analysis be performed between groups.

FoldEnrich is used for differential analysis of peaks in different experimental groups (the ratio of RPM value in group A to group B). Differerntial binding sites are analyzed in different groups by finding differential peak when the ratio of FoldEnrich is larger than 2. We can get the differential binding sites related genes to do the follow-up annotation and enrichment analysis.

8.1　Summary of differential peaks

The Differential Expression Venn diagram illustrates the differential peaks of the compared groups. The sum of all the numbers in the circle represents the total number in the compared groups, and the overlapping part is the common peak number between two compared groups., as shown in the figure below.

KD_vs_CTRL

Figure 8-1-1　Venn diagram of differential peak in comparable groups

Different colors represent different compared groups (A or B). The sum of the number in the pie is the peak number in this compared group, and the overlapping part is the number of common genes among comparison groups.

KD_vs_CTRL

Figure 8-1-2　Venn diagram of differential peak in comparable gruops

Since some compared groups (eg. A or B) may include multiple experimental groups (eg. A1, A2, A3), here show the number of peak among experimental groups within one compared group (eg. A1&A2, A1&A3). Different colors represent different experimental groups (A1 or A2). The sum of the number in the pie is the peak number of this experimental group. Single color is group specific peak number. Overlapping part is the common peak number between two experimental groups.

Novogene Co.,Ltd.

8.2　Enrichment level analysis among different samples

RPM value of the peak in different samples (the ratio of 1M reads that enriched to the peak in a single sample) is used to do clustering analysis to determine the enrichment pattern of the same peak in different samples or the enrichment difference of different peaks in the same sample. At the same time, Enrichment comparison between IP and Input in the group can show the peak enrichment in IP experiment. Hierarchical clustering is done based on RPM value in different samples. Different colors represent different clustering information. The more close between the IP sample enrichment, the more similar function or biological process they have.

KD_vs_CTRL_heatmap KD_vs_CTRL_corralation

Figure 8-2　Enrichment analysis among samples

Top is the clustering among different samples in comparable groups.

Bottom is the correlation among different samples in comparable groups.

8.3 Reads density distribution among comparable groups for the annotated peaks

RPM is the number of reads enriched in peak area per million reads. Boxplot is used to show reads enrichment in the annotated peak area in different samples.

KD_vs_CTRL

Figure 8-3 Boxplot of RPM value among samples in different groups

Novogene Co.,Ltd.

8.4　Enrichment level of different experimental peak analysis

FoldEnrich value (the Ratio of RPM value in IP to RPM in input) is used to do the cluster analysis of the peaks from different experimental group. It is used to determine the enrichment pattern of the same peak in different experimental groups, or the enrichment difference of different peaks in the same experimental group. At the same time, according to the comparison of peak enrichment among different experimental groups, we can show the different enrichment of peak between experimental groups, so as to find the different enrichment peak, which is the differential protein binding sites. Hierarchical clustering analysis is done based on the FoldEnrich value from different experimental groups. Different colors represent different clustering information.

KD_vs_CTRL_heatmap KD_vs_CTRL_corralation

Figure 8-4　Peak enrichment analysis among different experimental groups

Top is the clustering analysis of peaks from different experimental groups.

Bottom is the correlation among different experimental groups.

8.5　Reads density distribution among different experimental groups for the annotated peaks

Boxplot is used to show the different FoldEnrich (the ratio of RPM of IP to RPM of Input) i in the peak annotation region in different samples.

KD_vs_CTRL

Figure 8-5　Reads density distribution in peak annotation region from different experiment was shown by boxplot.

Novogene Co.,Ltd.

8.6　GO enrichment analysis of differential peak related genes

GO enrichment analysis for the differential peak related genes is shown below:

KD_vs_CTRL KD_vs_CTRL_mf KD_vs_CTRL_cc KD_vs_CTRL_bp KD_vs_CTRL_up KD_vs_CTRL_up_mf KD_vs_CTRL_up_cc KD_vs_CTRL_up_bp KD_vs_CTRL_down KD_vs_CTRL_down_mf KD_vs_CTRL_down_cc KD_vs_CTRL_down_bp

Figure 8-6　GO enrichment analysis for the differential peak genes

8.7　KEGG enrichment analysis of differential peak related genes

KEGG enrichment analysis for the differential peak related genes is shown below:

KD_vs_CTRL

Figure 8-7　KEGG enrichment analysis for differential peak genes

Novogene Co.,Ltd.

IV.Background

1　Experiment and library construction

ChIP-Seq is the efficient way to study protein-DNA interaction, transcription factor binding sites and histone modification in the whole genome (Landt and Marinov et al., 2012). The pipeline of ChIP-Seq and library construction is shown in the following figure. We provide services for every step started from library construction for the moment.

Normally, if you already know some target sequences, you can test the efficiency of the ChIP assay by qPCR before library construction. After sample quality control, library construction and quality control, we proceed it to the follow-up sequencing.

ChIP-seq library construction pipeline

Novogene Co.,Ltd.

2　Software list

Novogene Co.,Ltd.

3　References

Andrews S. (2010). FastQC: a quality control tool for high throughput sequence data. Available online at: http://www.bioinformatics.babraham.ac.uk/projects/fastqc

Ashburner, M. and C. A. Ball, et al. (2000). "Gene ontology: tool for the unification of biology. The Gene Ontology Consortium." Nat Genet 25 (1): 25-9.

Bailey, T. L. and N. Williams, et al. (2006). "MEME: discovering and analyzing DNA and protein sequence motifs." Nucleic Acids Res 34 (Web Server issue): W369-73.

Bailey, T. L. and M. Boden, et al. (2009). "MEME SUITE: tools for motif discovery and searching." Nucleic Acids Res 37 (Web Server issue): W202-8.

Bailey, T. and P. Krajewski, et al. (2013). "Practical guidelines for the comprehensive analysis of ChIP-seq data." PLoS Comput Biol 9 (11): e1003326.

Faust, G. G. and I. M. Hall (2014). "SAMBLASTER: fast duplicate marking and structural variant read extraction." Bioinformatics 30 (17): 2503-5.

Jiang, H. and R. Lei, et al. (2014). "Skewer: a fast and accurate adapter trimmer for next-generation sequencing paired-end reads." BMC Bioinformatics 15: 182.

Kanehisa, M. and S. Goto (2000). "KEGG: kyoto encyclopedia of genes and genomes." Nucleic Acids Res 28 (1): 27-30.

Kanehisa M., M. Araki, et al. (2008). KEGG for linking genomes to life and the environment. Nucleic acids research.(KEGG)

Kent, W. J. and A. S. Zweig, et al. (2010). "BigWig and BigBed: enabling browsing of large distributed datasets." Bioinformatics 26 (17): 2204-7. Available online at: http://hgdownload.cse.ucsc.edu/admin/exe/linux.x86_64/

Kharchenko, P. V. and M. Y. Tolstorukov, et al. (2008). "Design and analysis of ChIP-seq experiments for DNA-binding proteins." Nat Biotechnol 26 (12): 1351-9. Available online at: http://compbio.med.harvard.edu/Supplements/ChIP-seq/

Li, H. and R. Durbin (2009). "Fast and accurate short read alignment with Burrows-Wheeler transform." Bioinformatics 25 (14): 1754-60.

Li, H. and J. Ruan, et al. (2008). "Mapping short DNA sequencing reads and calling variants using mapping quality scores." Genome Res 18 (11): 1851-8.

Li, H. and B. Handsaker, et al. (2009). "The Sequence Alignment/Map format and SAMtools." Bioinformatics 25 (16): 2078-9.

Landt, S. G. and G. K. Marinov, et al. (2012). "ChIP-seq guidelines and practices of the ENCODE and modENCODE consortia." Genome Res 22 (9): 1813-31.

Mao, X. and T. Cai, et al. (2005). "Automated genome annotation and pathway identification using the KEGG Orthology (KO) as a controlled vocabulary." Bioinformatics 21 (19): 3787-93.

Nicol, J. W. and G. A. Helt, et al. (2009). "The Integrated Genome Browser: free software for distribution and exploration of genome-scale datasets." Bioinformatics 25 (20): 2730-1. Available online at: http://bioviz.org/igb/index.html

Peter J.Park(2009). ChIP-seq: advantages and challenges of a maturing technology. Nature Reviews Genetics 10, 669-679.

Ramirez, F. and F. Dundar, et al. (2014). "deepTools: a flexible platform for exploring deep-sequencing data." Nucleic Acids Res 42 (Web Server issue): W187-91.

R Core Team (2015). R: A Language and Environment for Statistical Computing. Available online at: https://www.r-project.org/

Shirley Pepke, Barbara Wold and Ali Mortazavi (2009).Computation for ChIP-seq and RNA-seq studies. Nature methods, VOL.6 NO.11s

Salmon-Divon, M. and H. Dvinge, et al. (2010). "PeakAnalyzer: genome-wide annotation of chromatin binding and modification loci." BMC Bioinformatics 11: 415.

Shobhit Gupta, JA Stamatoyannopolous, Timothy Bailey and William Stafford Noble, "Quantifying similarity between motifs", Genome Biology, 8(2):R24, 2007.

Yong Zhang,Tao Liu et al. (2008). Model-based Analysis of ChIP-Seq (MACS). Genome Biology, 9:R137

Young, M. D. and M. J. Wakefield, et al. (2010). "Gene ontology analysis for RNA-seq: accounting for selection bias." Genome Biol 11 (2): R14. Available online at: https://bioconductor.org/packages/release/bioc/html/goseq.html

Thorvaldsdottir, H. and J. T. Robinson, et al. (2013). "Integrative Genomics Viewer (IGV): high-performance genomics data visualization and exploration." Brief Bioinform 14 (2): 178-92. Available online at: https://www.broadinstitute.org/igv/

Timothy L. Bailey and Charles Elkan, "Fitting a mixture model by expectation maximization to discover motifs in biopolymers", Proceedings of the Second International Conference on Intelligent Systems for Molecular Biology, pp. 28-36, AAAI Press, Menlo Park, California, 1994.

Timothy L. Bailey, "DREME: Motif discovery in transcription factor ChIP-seq data", Bioinformatics, 27(12):1653-1659, 2011.

Young M D, Wakefield M J, Smyth G K, et al. (2010).Gene ontology analysis for RNA-seq: accounting for selection bias. Genome Biology, doi:10.1186/gb-2010-11-2-r14. (GOseq)

Kanehisa, M., M. Araki, et al. (2008). KEGG for linking genomes to life and the environment. Nucleic acids research.(KEGG)