libraries and construct a frequency distribution of the resulting
values (see Fig. 2). Note that after running the “collapsing_peaks” script, the number of 5
0 ends is greatly reduced and the
data can be manipulated with programs such as Excel and
Graphpad (see Note 17). We expect to see a bimodal distribution, in which the left peak represents cleavage sites and the
right peak represents TSSs. If the left peak is centered at zero,
the cleavage sites are equally represented in both library types
and no further normalization is necessary. If the left peak
deviates substantially from zero, choose a correction factor
and multiply all of the converted library normalization factors
by this value. Then, repeat the normalization in Subheading
3.10, step 4 using these new normalization factors. Visualize
the frequency distribution of log2 ratios again to determine if
the correction factor was appropriate and repeat this process as
needed until the left peak in the frequency distribution is
centered on or near zero.
3. Filter 5
0 ends based on their abundance relative to local expression library coverage (optional). Using strand-specific bam files
from an RNAseq expression library, Bedtools can report coverage at each coordinate with a command such as.
bedtools genomecov -d -pc -ibam Expression_library_sorted.
minus.bam -g NC_008596.chrom.sizes > Expression_library_sorted_minus_fragment.coverage.
With expression library coverage files for each strand in
hand, edit the script Mean_coverage_over_intervals.py as
needed and use it to determine the mean expression library
coverage in a desired window relative to each 5
0 end. We
recommend a 100 nt window upstream of each 5
0 end. Calculate the coverage ratio of each 5
0 end relative to the preceding
100 nt window in the expression library. You may then empirically choose a cutoff depending on your experimental goals (for
example, identification of major TSSs vs identification of sites
of low-level transcription initiation).
3.12 Distinguishing
TSSs from
Cleavage Sites
A frequency distribution of the log2 ratios of 5
0 end coverage in
converted/nonconverted libraries is expected to reveal a bimodal
distribution, in which a population of 5
0 ends with a log2 ratio
centered around zero represents cleavage sites and a population
with higher log2 ratios represents TSSs (see Fig. 2). We will use
the R package mixsmsn to perform mixture modeling to estimate
membership in each of these two populations. In the following
example, a text file contains a column titled “log2_phos_over_mock,” which reports the log2 ratios in converted/nonconverted
libraries for each 5
0 end:
528
M. Carla Martini et al.
values (see Fig. 2). Note that after running the “collapsing_peaks” script, the number of 5
0 ends is greatly reduced and the
data can be manipulated with programs such as Excel and
Graphpad (see Note 17). We expect to see a bimodal distribution, in which the left peak represents cleavage sites and the
right peak represents TSSs. If the left peak is centered at zero,
the cleavage sites are equally represented in both library types
and no further normalization is necessary. If the left peak
deviates substantially from zero, choose a correction factor
and multiply all of the converted library normalization factors
by this value. Then, repeat the normalization in Subheading
3.10, step 4 using these new normalization factors. Visualize
the frequency distribution of log2 ratios again to determine if
the correction factor was appropriate and repeat this process as
needed until the left peak in the frequency distribution is
centered on or near zero.
3. Filter 5
0 ends based on their abundance relative to local expression library coverage (optional). Using strand-specific bam files
from an RNAseq expression library, Bedtools can report coverage at each coordinate with a command such as.
bedtools genomecov -d -pc -ibam Expression_library_sorted.
minus.bam -g NC_008596.chrom.sizes > Expression_library_sorted_minus_fragment.coverage.
With expression library coverage files for each strand in
hand, edit the script Mean_coverage_over_intervals.py as
needed and use it to determine the mean expression library
coverage in a desired window relative to each 5
0 end. We
recommend a 100 nt window upstream of each 5
0 end. Calculate the coverage ratio of each 5
0 end relative to the preceding
100 nt window in the expression library. You may then empirically choose a cutoff depending on your experimental goals (for
example, identification of major TSSs vs identification of sites
of low-level transcription initiation).
3.12 Distinguishing
TSSs from
Cleavage Sites
A frequency distribution of the log2 ratios of 5
0 end coverage in
converted/nonconverted libraries is expected to reveal a bimodal
distribution, in which a population of 5
0 ends with a log2 ratio
centered around zero represents cleavage sites and a population
with higher log2 ratios represents TSSs (see Fig. 2). We will use
the R package mixsmsn to perform mixture modeling to estimate
membership in each of these two populations. In the following
example, a text file contains a column titled “log2_phos_over_mock,” which reports the log2 ratios in converted/nonconverted
libraries for each 5
0 end:
528
M. Carla Martini et al.