rna

CLASS1  = ["BMAL1-ZT10_rep1", "BMAL1-ZT10_rep2"]
CLASS1_PATHS = ["/mnt/Storage/data/exam_data/BMAL1/BMAL1-ZT10_rep1.bam", "/mnt/Storage/data/exam_data/BMAL1/BMAL1-ZT10_rep2.bam"]

CLASS2  = ["BMAL1-ZT22_rep1"]
CLASS2_PATHS = ["/mnt/Storage/data/exam_data/BMAL1/BMAL1-ZT22_rep1.bam"]

GTF = "/mnt/Storage/data/RefGene/mm9.gtf"


SAMPLES = CLASS1 + CLASS2
SAMPLE_PATHS = CLASS1_PATHS + CLASS2_PATHS

rule All:
     input: "diffexp/{A}#vs#{B}".format(A=",".join(CLASS1), B=",".join(CLASS2)), \
            expand("assembly/{sample}/transcripts.gtf", dir="assembly/{sample}", sample=SAMPLES)

rule SortBamFile:
     input: SAMPLE_PATHS
     output: protected("mapped/{sample}.bam".format(sample=sample)) for sample in SAMPLES
     message: "Sorting original BAM files.."
     log: "logs/sort.log"
     run:
     for one_input,one_output in zip(input, output):
          shell("samtools sort -m 5000000000 {0} {1}".format(one_input, one_output[:-4]))

rule CufflinksAssembly:
     input: "mapped/{sample}.bam"
     output: "assembly/{sample}/transcripts.gtf", dir="assembly/{sample}"
     log: "logs/assembly.log"
     threads: 4
     shell: "cufflinks --num-threads {threads} -o {output.dir} --frag-bias-correct {GTF} {input} -q"

rule GfoldCount:
     input: "mapped/{sample}.bam"
     output: "diffexp/{sample}.cnt"
     log: "logs/diffexp.log"
     shell: "samtools view {input} | gfold count -ann {GTF} -annf GTF -tag stdin -o {output}"

rule GfoldDiff:
     input: expand("diffexp/{sample}.cnt", sample=SAMPLES)
     output: "diffexp/{A}#vs#{B}".format(A=",".join(CLASS1), B=",".join(CLASS2))
     log: "logs/diffexp.log"
     run:
     shell("gfold diff -s1 {cl1} -s2 {cl2} -suf .cnt -o {output}".format(
     cl1 = ",".join(expand("diffexp/{sample}",sample=CLASS1)),
     cl2 = ",".join(expand("diffexp/{sample}",sample=CLASS2)),
     output = output))