RnaSnake/Includes/ExprAnalysis.snake

TRINITY_HOME = config["TRINITY_HOME"]
SAMPLES_FILE = config["SAMPLES_FILE"]
ExprAnalysis_Dir = config["ExprAnalysis_Dir"]
pc_num = config["pc_num"]
de_method = config["de_method"]

localrules: enrich

rule expr_stat:
    input:
        config["Quantification_Dir"] + "/my.gene.TMM.EXPR.matrix"
    output:
        config["Quantification_Dir"] + "/my.gene.TMM.EXPR.matrix.boxplot.pdf",
        config["Quantification_Dir"] + "/my.gene.TMM.EXPR.matrix.density.pdf",
        config["Quantification_Dir"] + "/my.gene.TMM.EXPR.matrix.stat.txt"
    shell:
        "Rscript script/expressionStats.R -d {input} -o {input}" 

rule sample_cor:
    input:
        config["Quantification_Dir"] + "/my.gene.TMM.EXPR.matrix"
    output:
        ExprAnalysis_Dir + "/sample_cor/my.minRow10.sample_cor.dat",
        ExprAnalysis_Dir + "/sample_cor/my.minRow10.sample_cor_matrix.pdf"
    shell:
        "mkdir -p {ExprAnalysis_Dir}/sample_cor;"
        "{TRINITY_HOME}/Analysis/DifferentialExpression/PtR"
        " --matrix {input}"
        " --samples {SAMPLES_FILE}"
        " --min_rowSums 10"
        " --sample_cor_matrix"
        " --output {ExprAnalysis_Dir}/sample_cor/my"

rule pca:
    input:
        config["Quantification_Dir"] + "/my.gene.TMM.EXPR.matrix"
    output:
        ExprAnalysis_Dir + "/pca/my.minRow10.PCA.prcomp.scores",
        ExprAnalysis_Dir + "/pca/my.minRow10.prcomp.principal_components.pdf"
    shell:
        "mkdir -p {ExprAnalysis_Dir}/pca;"
        "{TRINITY_HOME}/Analysis/DifferentialExpression/PtR"
        " --matrix {input}"
        " --samples {SAMPLES_FILE}"
        " --min_rowSums 10"
        " --prin_comp {pc_num}"
        " --add_prin_comp_heatmaps 10"
        " --output {ExprAnalysis_Dir}/pca/my"

contrasts = config["contrasts"]
f = open(contrasts, 'r')
CONTRASTS = [i.strip().replace('\t', '_vs_') for i in f.readlines()]

rule deg:
    input:
        config["Quantification_Dir"] + "/my.gene.counts.matrix"
    output:
        expand(ExprAnalysis_Dir + "/deg/my.gene.counts.matrix.{vs}.DESeq2.DE_results", vs = CONTRASTS)
    shell:
        "perl {TRINITY_HOME}/Analysis/DifferentialExpression/run_DE_analysis.pl"
        " --matrix {input} --method {de_method}"
        " --samples_file {SAMPLES_FILE}"
        " --contrasts  {contrasts}"
        " --output {ExprAnalysis_Dir}/deg"

de_log2FoldChange = config["de_log2FoldChange"]
de_padj = config["de_padj"]
ExprAnalysis_Dir = config["ExprAnalysis_Dir"]

rule enrich:
    input:
        config["Annotation_Dir"] + "/anno_stat.txt",
        des = expand(ExprAnalysis_Dir + "/deg/my.gene.counts.matrix.{vs}.DESeq2.DE_results", vs = CONTRASTS)
    output:
        expand(ExprAnalysis_Dir + "/deg/my.gene.counts.matrix.{vs}.DESeq2.DE_results.ekp_results.txt", vs = CONTRASTS),
        expand(ExprAnalysis_Dir + "/deg/my.gene.counts.matrix.{vs}.DESeq2.DE_results.ego_results.txt", vs = CONTRASTS)
    run:
        for i in input.des:
            shell("cd " + config["Annotation_Dir"] + "; Rscript ../tools/emcp/enrich.R -d ../{i} --de_log2FoldChange {de_log2FoldChange} --de_padj {de_padj}; cd ..")