add Pearson Correlation, make plots pdf

workflow/main.nf

@@ -198,7 +198,7 @@ process experimentQC {
 
   output:
 
-  file '*.{png,npz}' into deepToolsStats
+  file '*.{pdf,npz}' into deepToolsStats
 
   script:
 

workflow/scripts/experiment_qc.py

@@ -76,10 +76,10 @@ def generate_read_summary(design, extension):
     return mbs_filename
 
 
-def check_correlation(mbs):
+def check_spearman_correlation(mbs):
     '''Check Spearman pairwise correlation of samples based on read counts.'''
 
-    spearman_filename = 'heatmap_SpearmanCorr.png'
+    spearman_filename = 'heatmap_SpearmanCorr.pdf'
     spearman_params = "--corMethod spearman --skipZero" + \
                     " --plotTitle \"Spearman Correlation of Read Counts\"" + \
                     " --whatToPlot heatmap --colorMap RdYlBu --plotNumbers"
@@ -95,6 +95,25 @@ def check_correlation(mbs):
     out, err = spearman_correlation.communicate()
 
 
+def check_pearson_correlation(mbs):
+    '''Check Pearson pairwise correlation of samples based on read counts.'''
+
+    pearson_filename = 'heatmap_PearsonCorr.pdf'
+    pearson_params = "--corMethod pearson --skipZero" + \
+                    " --plotTitle \"Pearson Correlation of Read Counts\"" + \
+                    " --whatToPlot heatmap --colorMap RdYlBu --plotNumbers"
+
+    pearson_command = (
+        "plotCorrelation -in %s %s -o %s"
+        % (mbs, pearson_params, pearson_filename)
+    )
+
+    logger.info("Running deeptools with %s", pearson_command)
+
+    pearson_correlation = subprocess.Popen(pearson_command, shell=True)
+    out, err = pearson_correlation.communicate()
+
+
 def check_coverage(design, extension):
     '''Asses the sequencing depth of samples.'''
 
@@ -166,7 +185,9 @@ def main():
 
     # Run correlation
     mbs_filename = generate_read_summary(design_df, extension)
-    check_correlation(mbs_filename)
+    check_spearman_correlation(mbs_filename)
+    check_pearson_correlation(mbs_filename)
+
 
     # Run coverage
     check_coverage(design_df, extension)

workflow/tests/test_experiment_qc.py

@@ -33,7 +33,8 @@ def test_check_update_controls(design_bam):
 @pytest.mark.singleend
 def test_experiment_qc_singleend():
     assert os.path.exists(os.path.join(test_output_path, 'sample_mbs.npz'))
-    assert os.path.exists(os.path.join(test_output_path, 'heatmap_SpearmanCorr.png'))
+    assert os.path.exists(os.path.join(test_output_path, 'heatmap_SpearmanCorr.pdf'))
+    assert os.path.exists(os.path.join(test_output_path, 'heatmap_PearsonCorr.pdf'))
     assert os.path.exists(os.path.join(test_output_path, 'coverage.png'))
     assert os.path.exists(os.path.join(test_output_path, 'ENCLB144FDT_fingerprint.png'))
     assert os.path.exists(os.path.join(test_output_path, 'ENCLB831RUI_fingerprint.png'))
@@ -41,7 +42,8 @@ def test_experiment_qc_singleend():
 @pytest.mark.pairdend
 def test_experiment_qc_pairedend():
     assert os.path.exists(os.path.join(test_output_path, 'sample_mbs.npz'))
-    assert os.path.exists(os.path.join(test_output_path, 'heatmap_SpearmanCorr.png'))
+    assert os.path.exists(os.path.join(test_output_path, 'heatmap_SpearmanCorr.pdf'))
+    assert os.path.exists(os.path.join(test_output_path, 'heatmap_PearsonCorr.pdf'))
     assert os.path.exists(os.path.join(test_output_path, 'coverage.png'))
     assert os.path.exists(os.path.join(test_output_path, 'ENCLB568IYX_fingerprint.png'))
     assert os.path.exists(os.path.join(test_output_path, 'ENCLB637LZP_fingerprint.png'))