2
$\begingroup$

I have started analyzing methylation (EM-seq) data for the first time ~0.8M positions. In this data set I have 27 samples of 20 patients. I want to perform a PCA of the dataset to check for possible batch effects or problems.

Some publications and guides use the raw M-values to plot it without taking into account the coverage of each position. If I do that, I can detect a previously unknown source of variation, which is already great.

Usually there is also the information of the coverage of each position (not the number of reads per sample). Example of the data, which some packages store in the same object:

Cov:
      A      B      C      D      E
     19     10     18     14     20
     13      8     20     12     12
     12      6      9     16      7
     10      7      9     16      8
      8      5      6     12      8

M-values:
      A      B      C      D      E
   9.88   0.00   0.00   0.00   0.00
   7.93   0.00  10.00   3.00   0.00
   9.00   1.98   5.94  13.92   4.97
   9.00   2.94   3.96   8.96   4.96
   2.96   4.00   1.98   4.92   4.96

As you can see the coverage of each position is not the same for each sample. However, I haven't found any specific method for calculating a PCA taking into account the coverage: there doesn't seem to be a method that normalizes the M-values by the coverage of each position.

I checked the distribution of coverage and 91% is covered by less than 14 reads, but some positions are covered by more than 300 reads.

What is the best way to normalize the position by coverage and plot the PCA ?

$\endgroup$
5
  • 1
    $\begingroup$ This is because of the statistical issue of Eigen value calculations - transformation are generally not welcome. In the changing age of unsupervised learning the concept of PCA is changing admittedly. You can superimpose coverage post PCA to assess with it is affecting the output, i.e. if it starts forming clusters there's a problem and your approach is justified. $\endgroup$
    – M__
    Commented Jan 26, 2023 at 16:31
  • $\begingroup$ I have experience in the analysis of EM-seq sequencing datasets but I honestly never really used M-values instead of classic fractions of methylation. One common approach is too remove highly covered positions (outliers) as it is often the results of sequencing or mapping errors. If you have a average coverage of ~ 14X, I think a threshold at a 100X is good. On the detection of methylation variation accross sample conditions, you won't need to deal with coverage variations a thes DSS soft (for DMC and DMR detection) works with raw counts. Not sure if it helps! $\endgroup$ Commented Feb 2, 2023 at 10:52
  • $\begingroup$ @PaulEndymion I also used the beta-values for a PCA (it revealed two other sample outliers). Interesting point about the mapping and the errors, I haven't actually mapped it so I don't know how well it was done. Are you suggesting that I don't need to perform PCA because they are already taken into account? I will look if removing those highly covered positions changes anything. Thanks! $\endgroup$
    – llrs
    Commented Feb 3, 2023 at 10:32
  • $\begingroup$ Yes mapping errors on high repetition/GC% regions is a common bias, especially with short reads. These regions have the potential to attract reads with sequencing errors. We would do a PCA after the preprocessing of raw counts from "stranded_CpG_reports" files. Preprocessing included normalization of counts and removing of known SNP falling on "G" or "C" positions in "CG" motifs as well as higly covered positions. However as @M__ said if you found clusters of samples with PCA you can't explain with your samples experiment conditions, you might have a problem with these samples. $\endgroup$ Commented Feb 3, 2023 at 12:51
  • 1
    $\begingroup$ Just to be clear on the PCA. You can perform PCA but more as part of an exploratory/descriptive analysis, I don't remember reading about differential methylation studies using only PCA for DMR or DMC detection (even though DNA methylation studies are tedious so there is room for new approaches). $\endgroup$ Commented Feb 3, 2023 at 13:06

0

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Browse other questions tagged or ask your own question.