From c06b562c43990272940fabb9ed9b523a9b5ddde0 Mon Sep 17 00:00:00 2001
From: mariabernard <maria.bernard@jouy.inra.fr>
Date: Wed, 19 Jun 2019 12:56:34 +0200
Subject: [PATCH] 1000RNASEQ chicken ASE: add bai inputs, correct phaser.py
commande line
---
.../ASE/rules/ASE_counter.smk | 21 ++++++++++++++++---
1 file changed, 18 insertions(+), 3 deletions(-)
diff --git a/Snakemake/1000RNASeq_chicken/ASE/rules/ASE_counter.smk b/Snakemake/1000RNASeq_chicken/ASE/rules/ASE_counter.smk
index fc7ac09..95cfa1f 100644
--- a/Snakemake/1000RNASeq_chicken/ASE/rules/ASE_counter.smk
+++ b/Snakemake/1000RNASeq_chicken/ASE/rules/ASE_counter.smk
@@ -22,11 +22,25 @@ def get_splitNCigar_bam(wildcards):
else :
raise Exception(wildcards.sample + " is sequenced in pair end & single end mode\n This workflow do not accept this case!\n")
+def get_splitNCigar_bai(wildcards):
+
+ sample_table = table[table["sample_name"] == wildcards.sample]
+
+ # paired end, filter on properly paired and rmdup
+ if len(sample_table["forward_read"]) == len(sample_table["reverse_read"].dropna()) :
+ return "Results/STAR_Aln_2/" + wildcards.sample + "_rg_genomic_pp_rmdup_uniq_split.bai"
+ # single end, no filter
+ elif len(sample_table["reverse_read"].dropna() ) == 0:
+ return "Results/STAR_Aln_2/" + wildcards.sample + "_rg_genomic_rmdup_uniq_split.bai"
+ else :
+ raise Exception(wildcards.sample + " is sequenced in pair end & single end mode\n This workflow do not accept this case!\n")
+
rule ASEReadCounter:
input:
ref = config["fasta_ref"],
vcf = 'Results/hard_filtering/' + os.path.basename(config['vcf']).replace('.vcf.gz','') + '_SNP_filtered.vcf.gz',
- bam = get_splitNCigar_bam
+ bam = get_splitNCigar_bam ,
+ bai = get_splitNCigar_bai
output:
cvs = 'Results/Counter/{sample}_ASEReadCounter.csv'
log:
@@ -44,7 +58,8 @@ rule ASEReadCounter:
rule phASER:
input:
vcf = 'Results/hard_filtering/' + os.path.basename(config['vcf']).replace('.vcf.gz','') + '_SNP_filtered.vcf.gz',
- bam = get_splitNCigar_bam
+ bam = get_splitNCigar_bam,
+ bai = get_splitNCigar_bai
output:
cvs = 'Results/Counter/{sample}_phASER.csv'
log:
@@ -55,5 +70,5 @@ rule phASER:
paired_option = lambda wildcards : '1' if len(table[table['sample_name'] == wildcards.sample]['forward_read']) == len(table[table['sample_name'] == wildcards.sample]['reverse_read'].dropna()) else '0'
shell:
"""
- python2.7 phaser.py --vcf {input.vcf} --bam {input.bam} --sample {wildcards.sample} --threads {threads} --o Results/Counter/{wildcards.sample}_phASER --paired_end {params.paired_option} --baseq {params.baseQuality} --mapq {params.mappingQuality}
+ python2.7 {config[bin_dir]}/phaser.py --vcf {input.vcf} --bam {input.bam} --sample {wildcards.sample} --threads {threads} --o Results/Counter/{wildcards.sample}_phASER --paired_end {params.paired_option} --baseq {params.baseQuality} --mapq {params.mappingQuality}
"""
--
GitLab