Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
Original file line number Diff line number Diff line change
Expand Up @@ -1172,9 +1172,11 @@
"output: out\n",
"task: trunk_workers = 1, trunk_size = job_size, walltime = walltime, mem = mem, cores = numThreads, tags = f'{step_name}_{_output:bn}'\n",
"bash: expand= \"${ }\", stderr = f'{_output:n}.stderr', stdout = f'{_output:n}.stdout'\n",
" # The paths go in a file: a genome-wide region list is too long for one command-line argument.\n",
" printf '%s\\n' ${_input:r} > \"${_output:n}.data_files.txt\"\n",
" Rscript ${modular_script_dir}/data_preprocessing/genotype/genotype_formatting.R \\\n",
" --step write_data_list \\\n",
" --data-files \"${\"::\".join([str(x) for x in _input])}\" \\\n",
" --data-files-list \"${_output:n}.data_files.txt\" \\\n",
" --ext \"${ext}\" \\\n",
" --output \"${_output}\" \\\n",
" --numThreads ${numThreads}"
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -611,9 +611,11 @@
"input: group_by = \"all\"\n",
"output: f'{cwd}/{name}.{step_name[:-2]}_files.txt'\n",
"bash: expand= \"${ }\", stderr = f'{_output:n}.stderr', stdout = f'{_output:n}.stdout', container = container, entrypoint = entrypoint\n",
" # The paths go in a file: a genome-wide region list is too long for the command line.\n",
" printf '%s\\n' ${_input:r} > \"${_output:n}.inputs.txt\"\n",
" Rscript ${modular_script_dir}/data_preprocessing/phenotype/phenotype_formatting.R --step phenotype_by_region_2 \\\n",
" --cwd \"${cwd}\" \\\n",
" --inputs ${_input} \\\n",
" --inputs-list \"${_output:n}.inputs.txt\" \\\n",
" --output \"${_output}\""
]
},
Expand Down
4 changes: 3 additions & 1 deletion code/SoS/mnm_analysis/mnm_methods/colocboost.ipynb
Original file line number Diff line number Diff line change
Expand Up @@ -497,7 +497,9 @@
"parameter: ld_meta_data = path()\n",
"# Optional region list whose LAST column lists region/gene IDs to analyze.\n",
"parameter: region_list = path()\n",
"# Region/gene IDs to analyze (space separated on the CLI).\n",
"# Region/gene IDs to analyze (space separated on the CLI). For thousands of IDs\n",
"# use --region-list instead: the IDs reach colocboost_manifest.R as one\n",
"# comma-joined argument, which Linux caps at 128 KiB.\n",
"parameter: region_name = []\n",
"# Optional per-gene association windows (BED: #chr start end ID) overriding the\n",
"# default cis_window for the GWAS extraction window.\n",
Expand Down
132 changes: 91 additions & 41 deletions code/SoS/mnm_analysis/mnm_methods/mnm_regression.ipynb

Large diffs are not rendered by default.

34 changes: 22 additions & 12 deletions code/SoS/mnm_analysis/mnm_postprocessing.ipynb
Original file line number Diff line number Diff line change
Expand Up @@ -97,7 +97,7 @@
"\n",
"The remaining workflows take these additional flags:\n",
"\n",
"- `--tsv_path output_postproc/*.lbf.tsv` (per-region LBF tables to concatenate; `susie_tsv_collapse`)\n",
"- `--tsv_path 'output_postproc/*.lbf.tsv'` (per-region LBF tables to concatenate; `susie_tsv_collapse`. Quote the wildcard so SoS expands it: expanded by the shell, a genome-wide set of files is too long for the command line.)\n",
"- `--annot_tibble <annotation_tibble.tsv>` (annotation tibble used to label variants; `susie_pip_landscape_plot` and `fsusie_extract_effect`)\n",
"- `--trait_to_select 1` (how many traits the credible-set UpSet diagram shows; `susie_upsetR_cs_plot`)\n",
"- `--SNP_list <variant_list.txt>` (variants to annotate; `tmp_annotatation_of_snps`)\n",
Expand Down Expand Up @@ -264,7 +264,7 @@
"source": [
"sos run pipeline/mnm_postprocessing.ipynb susie_tsv_collapse \\\n",
" --cwd output_postproc --study <study_name> \\\n",
" --tsv_path output_postproc/*.lbf.tsv \\\n",
" --tsv_path 'output_postproc/*.lbf.tsv' \\\n",
" --container '' -j1"
]
},
Expand All @@ -276,7 +276,7 @@
},
"source": [
"### Diagnostic plots\n",
"`susie_pip_landscape_plot` draws the PIP landscape; `susie_upsetR_plot` and `susie_upsetR_cs_plot` draw UpSet diagrams of credible-set overlap. The landscape plot needs an annotation tibble; `--trait_to_select` picks how many traits the credible-set diagram shows."
"`susie_pip_landscape_plot` draws the PIP landscape; `susie_upsetR_plot` and `susie_upsetR_cs_plot` draw UpSet diagrams of credible-set overlap. The landscape plot needs an annotation tibble; `--trait_to_select` picks how many traits the credible-set diagram shows. Quote the `--rds_path` wildcard so SoS expands it: expanded by the shell, a genome-wide set of files is too long for the command line."
]
},
{
Expand All @@ -300,13 +300,13 @@
"source": [
"sos run pipeline/mnm_postprocessing.ipynb susie_pip_landscape_plot \\\n",
" --cwd output_postproc --study <study_name> \\\n",
" --rds_path <fine_mapping_dir>/*.fine_mapping.rds \\\n",
" --rds_path '<fine_mapping_dir>/*.fine_mapping.rds' \\\n",
" --annot_tibble <annotation_tibble.tsv> \\\n",
" --container '' -j1\n",
"\n",
"sos run pipeline/mnm_postprocessing.ipynb susie_upsetR_cs_plot \\\n",
" --cwd output_postproc --study <study_name> \\\n",
" --rds_path <fine_mapping_dir>/*.fine_mapping.rds \\\n",
" --rds_path '<fine_mapping_dir>/*.fine_mapping.rds' \\\n",
" --trait_to_select 1 \\\n",
" --container '' -j1"
]
Expand Down Expand Up @@ -1098,8 +1098,10 @@
"output: f'{cwd}/{name}.UpSetR.pdf'\n",
"task: trunk_workers = 1, trunk_size = job_size, walltime = walltime, mem = mem, cores = numThreads, tags = f'{step_name}_{_output:bn}'\n",
"bash: expand = '${ }', stdout = f\"{_output:n}.stdout\", stderr = f\"{_output:n}.stderr\", container = container, entrypoint = entrypoint\n",
" # The paths go in a file: a genome-wide list is too long for the command line.\n",
" printf '%s\\n' ${_input:r} > \"${_output:n}.input.txt\"\n",
" Rscript ${modular_script_dir}/pecotmr_integration/fine_mapping_upset.R \\\n",
" --input ${\" \".join([str(x) for x in _input])} \\\n",
" --input-list \"${_output:n}.input.txt\" \\\n",
" --output ${_output}\n"
]
},
Expand All @@ -1123,8 +1125,10 @@
"output: f'{cwd}/{name}.UpSetR_cs.pdf'\n",
"task: trunk_workers = 1, trunk_size = job_size, walltime = walltime, mem = mem, cores = numThreads, tags = f'{step_name}_{_output:bn}'\n",
"bash: expand = '${ }', stdout = f\"{_output:n}.stdout\", stderr = f\"{_output:n}.stderr\", container = container, entrypoint = entrypoint\n",
" # The paths go in a file: a genome-wide list is too long for the command line.\n",
" printf '%s\\n' ${_input:r} > \"${_output:n}.input.txt\"\n",
" Rscript ${modular_script_dir}/pecotmr_integration/fine_mapping_upset.R \\\n",
" --input ${\" \".join([str(x) for x in _input])} \\\n",
" --input-list \"${_output:n}.input.txt\" \\\n",
" --output ${_output}\n"
]
},
Expand Down Expand Up @@ -1505,8 +1509,10 @@
"task: trunk_workers = 1, walltime = '1h', trunk_size = 1, mem = '16G', cores = 1, tags = f'{_output:bn}'\n",
"bash: expand = \"${ }\", container = container, stderr = f'{_output:n}.stderr', stdout = f'{_output:n}.stdout'\n",
" if [ -e \"${_output:ad}/${name}_cache/\" ]; then\n",
" sed 's/^chr/#chr/' `ls ${_output:ad}/${name}_cache/*tsv |head -n1` | head -n 1 > ${_output:an}.temp\n",
" tail -n +2 -q ${_output:ad}/${name}_cache/*.tsv >> ${_output:an}.temp\n",
" # find + xargs, not a glob: a genome-wide run has one cache table per region, too many for one command line.\n",
" first_table=$(find \"${_output:ad}/${name}_cache/\" -maxdepth 1 -type f -name '*tsv' | sort | head -n 1)\n",
" sed 's/^chr/#chr/' \"$first_table\" | head -n 1 > ${_output:an}.temp\n",
" find \"${_output:ad}/${name}_cache/\" -maxdepth 1 -type f -name '*.tsv' -print0 | sort -z | xargs -0 -r tail -n +2 -q >> ${_output:an}.temp\n",
" error_files=$(find \"${_output:ad}/${name}_cache/\" -type f -name \"*_error\")\n",
"\n",
" if [[ -n $error_files ]]; then\n",
Expand Down Expand Up @@ -1625,8 +1631,10 @@
"parameter: remove_cache = False\n",
"output: f\"{cache_path:d}/{output_file}\" \n",
"bash: expand = \"${ }\", container = container, stderr = f'{_input:n}.stderr', stdout = f'{_input:n}.stdout' \n",
" head -n 1 -q ${cache_path}/*.tsv | sed 's/^chr/#chr/' | head -n 1 > ${_output}\n",
" tail -n +2 -q ${cache_path}/*.tsv >> ${_output}\n",
" # find + xargs, not a glob: one cache table per region is too many for one command line.\n",
" first_table=$(find \"${cache_path}\" -maxdepth 1 -type f -name '*.tsv' | sort | head -n 1)\n",
" head -n 1 \"$first_table\" | sed 's/^chr/#chr/' > ${_output}\n",
" find \"${cache_path}\" -maxdepth 1 -type f -name '*.tsv' -print0 | sort -z | xargs -0 -r tail -n +2 -q >> ${_output}\n",
" if ${\"true\" if remove_cache else \"false\"}; then\n",
" rm -rf ${cache_path}\n",
" fi"
Expand Down Expand Up @@ -1708,8 +1716,10 @@
"output: f\"{cwd}/summary/{name}.{qtl_type}{_tag}.top_loci.bed.gz\"\n",
"task: trunk_workers = 1, walltime = '2h', trunk_size = 1, mem = '8G', cores = 1, tags = f'{_output:bn}'\n",
"bash: expand = \"${ }\", container = container, stderr = f'{_output:n}.stderr', stdout = f'{_output:n}.stdout', entrypoint=entrypoint\n",
" # The paths go in a file: a genome-wide list is too long for the command line.\n",
" printf '%s\\n' ${_input:r} > \"${_output:nn}.input.txt\"\n",
" Rscript ${modular_script_dir}/pecotmr_integration/fine_mapping_bed_merge.R \\\n",
" --input ${_input:r} \\\n",
" --input-list \"${_output:nn}.input.txt\" \\\n",
" --output ${_output:r}"
]
},
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -554,7 +554,10 @@
"output: f'{cwd}/{leafcutter_version}/{dataset}_intron_usage_perind.junc'\n",
"bash: container=container, expand= \"${ }\", stderr = f'{_output[0]:n}.stderr', stdout = f'{_output[0]:n}.stderr' \n",
" cd ${cwd}\n",
" realpath ${junc_path}/*${file_suffix} > ${_output}.tmp\n",
" # find, not a glob: a cohort's junction files are too many arguments for\n",
" # realpath. Hidden files are skipped, as the glob skipped them.\n",
" find \"${junc_path}/\" -mindepth 1 -maxdepth 1 -name \"*${file_suffix}\" ! -name '.*' -print0 \\\n",
" | sort -z | xargs -0 -r realpath > ${_output}.tmp\n",
" grep -v \"all\" ${_output}.tmp > ${_output}\n",
" rm ${_output}.tmp"
]
Expand Down
16 changes: 12 additions & 4 deletions code/SoS/molecular_phenotypes/calling/RNA_calling.ipynb
Original file line number Diff line number Diff line change
Expand Up @@ -954,11 +954,13 @@
"input: group_by = \"all\"\n",
"output: f'{cwd}/{sample_list:n}.trimmed.txt'\n",
"bash: expand = \"${ }\", stderr = f'{_output:n}.stderr', stdout = f'{_output:n}.stdout'\n",
" # The paths go in a file: a large cohort's list is too long for the command line.\n",
" printf '%s\\n' ${_input:r} > \"${_output:n}.input.txt\"\n",
" Rscript ${modular_script_dir}/molecular_phenotypes/calling/RNA_calling.R \\\n",
" --step fastp_manifest \\\n",
" --sample-list \"${sample_list}\" \\\n",
" --is-paired-end ${is_paired_end} \\\n",
" --input ${' '.join([f'\"{x}\"' for x in _input])} \\\n",
" --input-list \"${_output:n}.input.txt\" \\\n",
" --output \"${_output}\"\n"
]
},
Expand Down Expand Up @@ -1255,10 +1257,12 @@
"output: f'{cwd}/{sample_list:bn}.bam_file_list'\n",
"task: trunk_workers = 1, trunk_size = 1, walltime = walltime, mem = mem, cores = 1\n",
"bash: expand= \"${ }\", stderr = f'{_output:n}.stderr', stdout = f'{_output:n}.stdout'\n",
" # The paths go in a file: a large cohort's list is too long for the command line.\n",
" printf '%s\\n' ${_input:r} > \"${_output}.input.txt\"\n",
" Rscript ${modular_script_dir}/molecular_phenotypes/calling/RNA_calling.R \\\n",
" --step star_align_3 \\\n",
" --output \"${_output}\" \\\n",
" --input ${\" \".join([str(x) for x in _input])} \\\n",
" --input-list \"${_output}.input.txt\" \\\n",
" --sample-id ${\" \".join([str(x) for x in sample_id])} \\\n",
" --strand ${\" \".join([str(x) for x in strand])} \\\n",
" --var-vcf-file \"${varVCFfile}\"\n"
Expand Down Expand Up @@ -1360,11 +1364,13 @@
" f'{cwd}/{bam_list:bn}.rnaseqc.metrics.tsv'\n",
"task: trunk_workers = 1, trunk_size = job_size, walltime = walltime, mem = mem, cores = numThreads\n",
"bash: expand= \"${ }\", stderr = f'{_output[0]:n}.stderr', stdout = f'{_output[0]:n}.stdout'\n",
" # The paths go in a file: a large cohort's list is too long for the command line.\n",
" printf '%s\\n' ${_input:r} > \"${cwd}/${bam_list:bn}.rnaseqc.input.txt\"\n",
" Rscript ${modular_script_dir}/molecular_phenotypes/calling/RNA_calling.R \\\n",
" --step rnaseqc_merge \\\n",
" --cwd \"${cwd}\" \\\n",
" --name \"${bam_list:bn}\" \\\n",
" --input ${\" \".join([str(x) for x in _input])}\n"
" --input-list \"${cwd}/${bam_list:bn}.rnaseqc.input.txt\"\n"
]
},
{
Expand Down Expand Up @@ -1457,11 +1463,13 @@
"\n",
"task: trunk_workers = 1, walltime = walltime, mem = mem, cores = numThreads, trunk_size = job_size\n",
"bash: expand= \"${ }\", stderr = f'{_output[0]:n}.stderr', stdout = f'{_output[0]:n}.stdout'\n",
" # The paths go in a file: a large cohort's list is too long for the command line.\n",
" printf '%s\\n' ${_input:r} > \"${cwd}/${bam_list:bn}.rsem.input.txt\"\n",
" Rscript ${modular_script_dir}/molecular_phenotypes/calling/RNA_calling.R \\\n",
" --step rsem_call_2 \\\n",
" --cwd \"${cwd}\" \\\n",
" --name \"${bam_list:bn}\" \\\n",
" --input ${\" \".join([str(x) for x in _input])}\n"
" --input-list \"${cwd}/${bam_list:bn}.rsem.input.txt\"\n"
]
},
{
Expand Down
12 changes: 9 additions & 3 deletions code/SoS/multivariate_genome/MASH/mash_posterior.ipynb
Original file line number Diff line number Diff line change
Expand Up @@ -764,8 +764,10 @@
"input: group_by = \"all\"\n",
"output: f\"{cwd}/{name}.posterior_sum.csv\"\n",
"bash: expand = \"${ }\", stderr = f'{_output:n}.stderr', stdout = f'{_output:n}.stdout', container = container\n",
" # The paths go in a file: a genome-wide list is too long for the command line.\n",
" printf '%s\\n' ${_input:r} > \"${_output:n}.contrast.txt\"\n",
" Rscript ${modular_script_dir}/pecotmr_integration/mash_posterior_contrast_summary.R \\\n",
" --contrast ${_input} \\\n",
" --contrast-list \"${_output:n}.contrast.txt\" \\\n",
" --cells ${','.join(cells)} \\\n",
" --p-cutoff ${p_cut} \\\n",
" --output ${_output}"
Expand Down Expand Up @@ -886,8 +888,10 @@
"input: group_by = \"all\"\n",
"output: f\"{cwd}/{name}.{step_name}.feature_score_sum.csv\"\n",
"bash: expand = \"${ }\", stderr = f'{_output:n}.stderr', stdout = f'{_output:n}.stdout', container = container\n",
" # The paths go in a file: a genome-wide list is too long for the command line.\n",
" printf '%s\\n' ${_input:r} > \"${_output:n}.scores.txt\"\n",
" Rscript ${modular_script_dir}/pecotmr_integration/mash_feature_score_merge.R \\\n",
" --scores ${_input} \\\n",
" --scores-list \"${_output:n}.scores.txt\" \\\n",
" --output ${_output}"
]
},
Expand Down Expand Up @@ -983,8 +987,10 @@
"input: group_by = \"all\"\n",
"output: f\"{cwd}/{name}.{step_name}.feature_score_sum.csv\"\n",
"bash: expand = \"${ }\", stderr = f'{_output:n}.stderr', stdout = f'{_output:n}.stdout', container = container\n",
" # The paths go in a file: a genome-wide list is too long for the command line.\n",
" printf '%s\\n' ${_input:r} > \"${_output:n}.scores.txt\"\n",
" Rscript ${modular_script_dir}/pecotmr_integration/mash_feature_score_merge.R \\\n",
" --scores ${_input} \\\n",
" --scores-list \"${_output:n}.scores.txt\" \\\n",
" --output ${_output}"
]
}
Expand Down
13 changes: 9 additions & 4 deletions code/SoS/multivariate_genome/MASH/mash_preprocessing.ipynb
Original file line number Diff line number Diff line change
Expand Up @@ -538,14 +538,17 @@
"with open(fine_mapping_meta) as _f:\n",
" _rows = list(csv.DictReader(_f, delimiter = '\\t'))\n",
"_fmr_paths = [r[finemapping_column] for r in _rows]\n",
"_region_ids = \",\".join(str(r[region_column]) for r in _rows)\n",
"_region_ids = [str(r[region_column]) for r in _rows]\n",
"input: _fmr_paths\n",
"output: f\"{cwd}/{name}.mash_input.rds\"\n",
"task: trunk_workers = 1, trunk_size = job_size, walltime = walltime, mem = mem, cores = numThreads, tags = f'{step_name}_{_output:bn}'\n",
"bash: expand = \"${ }\", stderr = f'{_output:n}.stderr', stdout = f'{_output:n}.stdout', container = container\n",
" # The paths and region ids go in files: genome-wide lists are too long for the command line.\n",
" printf '%s\\n' ${_input:r} > \"${_output:n}.objects.txt\"\n",
" printf '%s\\n' ${\" \".join(repr(x) for x in _region_ids)} > \"${_output:n}.region_ids.txt\"\n",
" Rscript ${modular_script_dir}/pecotmr_integration/mash_preprocessing.R \\\n",
" --objects ${_input} \\\n",
" --region-ids ${_region_ids} \\\n",
" --objects-list \"${_output:n}.objects.txt\" \\\n",
" --region-ids-list \"${_output:n}.region_ids.txt\" \\\n",
" --coverage ${coverage} \\\n",
" --n-random ${n_random} \\\n",
" --n-null ${n_null} \\\n",
Expand Down Expand Up @@ -622,8 +625,10 @@
"output: f\"{cwd}/{name}.mash_input.rds\"\n",
"task: trunk_workers = 1, trunk_size = job_size, walltime = walltime, mem = mem, cores = numThreads, tags = f'{_output:bn}'\n",
"bash: expand = \"${ }\", stderr = f'{_output:n}.stderr', stdout = f'{_output:n}.stdout', container = container\n",
" # The paths go in a file: a genome-wide list is too long for the command line.\n",
" printf '%s\\n' ${_input:r} > \"${_output:n}.objects.txt\"\n",
" Rscript ${modular_script_dir}/pecotmr_integration/mash_preprocessing.R \\\n",
" --objects ${_input} \\\n",
" --objects-list \"${_output:n}.objects.txt\" \\\n",
" --n-random ${n_random} \\\n",
" --n-null ${n_null} \\\n",
" --exclude-condition \"${','.join(exclude_condition)}\" \\\n",
Expand Down
9 changes: 7 additions & 2 deletions code/SoS/pecotmr_integration/SuSiE_enloc.ipynb
Original file line number Diff line number Diff line change
Expand Up @@ -435,6 +435,8 @@
"parameter: region_list = path()\n",
"# Optional: if a region name is provided \n",
"# the analysis would be focused on the union of provides region list and region names\n",
"# For thousands of regions use --region-list instead: the names reach\n",
"# enloc_manifest.R as one comma-joined argument, which Linux caps at 128 KiB.\n",
"parameter: region_name = []\n",
"# It is required to input the name of the analysis\n",
"parameter: name = f\"{xqtl_meta_data:bnn}.{gwas_meta_data:bnn}\"\n",
Expand Down Expand Up @@ -571,9 +573,12 @@
" # The xQTL and GWAS sides are already S4 FineMappingResult collections\n",
" # (fine_mapping output); qtl_enrichment.R combines its multi-file inputs\n",
" # itself, so the legacy legacy_enloc_finemap_convert.R step is gone.\n",
" # The paths go in files: a condition's genome-wide list is too long for the command line.\n",
" printf '%s\\n' ${_input:r} > \"${_output:n}.qtl_fine_mapping.txt\"\n",
" printf '%s\\n' ${\" \".join(repr(x) for x in _jobs[\"gwas_files\"].split(\",\"))} > \"${_output:n}.gwas_fine_mapping.txt\"\n",
" Rscript ${modular_script_dir}/pecotmr_integration/qtl_enrichment.R \\\n",
" --qtl-fine-mapping ${\" \".join([str(x) for x in _input])} \\\n",
" --gwas-fine-mapping ${_jobs[\"gwas_files\"].replace(\",\", \" \")} \\\n",
" --qtl-fine-mapping-list \"${_output:n}.qtl_fine_mapping.txt\" \\\n",
" --gwas-fine-mapping-list \"${_output:n}.gwas_fine_mapping.txt\" \\\n",
" --ncore ${numThreads} \\\n",
" ${seed_arg} \\\n",
" --output ${_output}\n"
Expand Down
Loading
Loading