From 9a450c94bd9203056bafe1385bb870189d125e3e Mon Sep 17 00:00:00 2001 From: Floreal Cabanettes <floreal.cabanettes@inra.fr> Date: Wed, 6 Jun 2018 17:52:21 +0200 Subject: [PATCH] Add align after build_pop --- cnvpipelines.py | 22 +++++++++++++++++++--- snakecnv/align.snk | 6 +++--- snakecnv/popsim.snk | 12 +++++++----- 3 files changed, 29 insertions(+), 11 deletions(-) diff --git a/cnvpipelines.py b/cnvpipelines.py index da873d2..b0e8fc5 100755 --- a/cnvpipelines.py +++ b/cnvpipelines.py @@ -705,13 +705,13 @@ class CnvPipeline: raise FileExistsError(file + "already exists but is not a link") os.symlink(file, final_link) - final_sample_file = os.path.join(self.wdir, "samples.list") + final_sample_file = os.path.join(self.wdir, "samples.yml") self.parse_align_sample_file(samples, final_sample_file) config = { "wdir": self.wdir, "reference": os.path.join(ref_dir, os.path.basename(reference)), - "sample_file": final_sample_file + "sample_file_align": final_sample_file } with open(self.config_file, "w") as config_file: @@ -753,6 +753,17 @@ class CnvPipeline: if not os.path.exists(dest): os.symlink(source, dest) + @staticmethod + def _create_samples_fastq_simulaton(nb_inds, smple_file): + samples = {} + for i in range(1, nb_inds+1): + samples["indiv%d" % i] = { + "reads1": os.path.join("pop", "indiv%d_1.fq.gz" % i), + "reads2": os.path.join("pop", "indiv%d_2.fq.gz" % i) + } + with open(smple_file, "w") as samples_f: + yaml.dump(samples, samples_f, default_flow_style=False) + def run_simulation(self, nb_inds, reference, nstretches, sv_list, coverage, force_polymorphism, haploid, proba_del, proba_inv, read_len, insert_len_mean, insert_len_sd, min_deletions, min_inversions, max_try, genotypes, force_wdir=False, **kwargs): @@ -799,6 +810,10 @@ class CnvPipeline: ref_final = os.path.join(self.wdir, os.path.basename(reference)) self._link(reference, ref_final) + # Create sample file for fastq files: + samples_file_fq = os.path.join(self.wdir, "samples.yml") + self._create_samples_fastq_simulaton(nb_inds, samples_file_fq) + config = { "wdir": self.wdir, "nb_inds": nb_inds, @@ -813,7 +828,8 @@ class CnvPipeline: "insert_len_sd": insert_len_sd, "min_deletions": min_deletions, "min_inversions": min_inversions, - "max_try": max_try + "max_try": max_try, + "sample_file_align": samples_file_fq } nstretches_final = None diff --git a/snakecnv/align.snk b/snakecnv/align.snk index c051364..eaba3ab 100644 --- a/snakecnv/align.snk +++ b/snakecnv/align.snk @@ -13,7 +13,7 @@ from align_snakemake_utils import SnakemakeUtils # Get from config: -samples = get_samples_align(config['sample_file']) +samples = get_samples_align(config['sample_file_align']) reference = config['reference'] wdir = config['wdir'] @@ -23,10 +23,10 @@ snakemake_utils = SnakemakeUtils(samples) include: "tools/threads.snk" -localrules: all +localrules: all_align -rule all: +rule all_align: input: expand("bams/{sample}.bam", sample=samples.keys()), expand("bams/{sample}.bam.bai", sample=samples.keys()) diff --git a/snakecnv/popsim.snk b/snakecnv/popsim.snk index 7458e7d..3a41aef 100644 --- a/snakecnv/popsim.snk +++ b/snakecnv/popsim.snk @@ -54,7 +54,7 @@ def buildpop_inputs(): def buildpop_outputs(): outputs = [] for n in range(1, NB_INDS+1): - outputs += [os.path.join("pop", "indiv%d_1.fq.gz" % n), os.path.join("pop", "indiv%d_2.fq.gz" % n)] + outputs += [os.path.abspath(os.path.join("pop", "indiv%d_1.fq.gz" % n)), os.path.abspath(os.path.join("pop", "indiv%d_2.fq.gz" % n))] return outputs bp_inputs, bp_params = buildpop_inputs() @@ -62,12 +62,14 @@ bp_outputs = buildpop_outputs() include: "tools/threads.snk" +include: "align.snk" -localrules: all +localrules: all_popsim -rule all: +rule all_popsim: input: - *bp_outputs + expand("bams/{sample}.bam", sample=["indiv%d" % i for i in range(1, NB_INDS+1)]) + #*bp_outputs rule buildpop: input: @@ -118,4 +120,4 @@ rule buildpop: command.append("-f") if params.haploid: command.append("-a") - shell(" ".join(command) + " 1> %s 2> %s" % (log.stdout, log.stderr)) \ No newline at end of file + shell(" ".join(command) + " 1> %s 2> %s" % (log.stdout, log.stderr)) -- GitLab