diff --git a/.github/pull_request_template.md b/.github/pull_request_template.md new file mode 100644 index 0000000..16d18a6 --- /dev/null +++ b/.github/pull_request_template.md @@ -0,0 +1,7 @@ +## Pull request (PR) description: +> Provide a brief description of the PR here. +> +> It does not need to list what each individual commit does, but rather provide a high-level overview of the PR. + +## Issues resolved: +> Any associated issues should be provided here using the format: 'Resolves #' (e.g. Resolves #1). If there are none, use 'n/a' instead. diff --git a/.github/workflows/pr-ci-report.yml b/.github/workflows/pr-ci-report.yml new file mode 100644 index 0000000..1ed6fe6 --- /dev/null +++ b/.github/workflows/pr-ci-report.yml @@ -0,0 +1,249 @@ +name: Pull Request CI Report + +on: + pull_request: + +permissions: + contents: read + pull-requests: write + +jobs: + python-tests: + name: Python ${{ matrix.python-version }} (${{ matrix.os }}) + runs-on: ${{ matrix.os }} + if: github.event.pull_request.head.repo.full_name == github.repository + strategy: + fail-fast: false + matrix: + os: [ubuntu-latest, macos-latest] + python-version: ["3.14", "3.15"] + steps: + - uses: actions/checkout@v7 + + - name: OS version (Linux) + if: runner.os == 'Linux' + shell: bash + run: | + . /etc/os-release + echo "OS_NAME=$NAME" >> "$GITHUB_ENV" + echo "OS_VER=$VERSION_ID" >> "$GITHUB_ENV" + + - name: OS version (macOS) + if: runner.os == 'macOS' + shell: bash + run: | + echo "OS_NAME=macOS" >> "$GITHUB_ENV" + echo "OS_VER=$(sw_vers -productVersion)" >> "$GITHUB_ENV" + + - name: Install Linux Qt runtime dependencies + if: runner.os == 'Linux' + shell: bash + run: | + sudo apt-get update + sudo apt-get install -y \ + libegl1 \ + libgl1 \ + libdbus-1-3 \ + libxkbcommon-x11-0 \ + libxcb-cursor0 \ + libxcb-icccm4 \ + libxcb-image0 \ + libxcb-keysyms1 \ + libxcb-randr0 \ + libxcb-render-util0 \ + libxcb-shape0 \ + libxcb-xinerama0 \ + libxcb-xfixes0 \ + xvfb + + - name: Install uv + uses: astral-sh/setup-uv@ae62891fec2bb8e7d6c99fc78c9fec3a63790f8d + with: + python-version: ${{ matrix.python-version }} + + - name: Sync dependencies + run: uv sync --all-extras + + - name: Run pytest + id: pytest + continue-on-error: true + run: uv run pytest + + - name: Build package + id: build + continue-on-error: true + run: uv build + + - name: Add venv bin to PATH + shell: bash + run: echo "$PWD/.venv/bin" >> "$GITHUB_PATH" + + - name: Check CLI entry point + id: cli + continue-on-error: true + run: comms --help + + - name: Write result file + if: always() + shell: bash + run: | + mkdir -p result + cat > result/result.json <> "$GITHUB_ENV" + echo "OS_VER=$VERSION_ID" >> "$GITHUB_ENV" + echo "USE_BUNDLED_LIBUV=1" >> "$GITHUB_ENV" + + - name: OS version (macOS) + if: runner.os == 'macOS' + shell: bash + run: | + echo "OS_NAME=macOS" >> "$GITHUB_ENV" + echo "OS_VER=$(sw_vers -productVersion)" >> "$GITHUB_ENV" + + - name: Install Linux runtime dependencies + if: runner.os == 'Linux' + shell: bash + run: | + sudo apt-get update + sudo apt-get install -y \ + libcurl4-openssl-dev \ + libfontconfig1-dev \ + libfreetype6-dev \ + libfribidi-dev \ + libharfbuzz-dev \ + libjpeg-dev \ + libpng-dev \ + libtiff5-dev \ + libwebp-dev + + - name: Set up R + uses: r-lib/actions/setup-r@v2 + with: + r-version: ${{ matrix.r-version }} + use-public-rspm: true + + - name: Install R dependencies + run: | + Rscript -e 'install.packages("testthat", repos = "https://cloud.r-project.org")' + Rscript src/comms/r/deps/install_deps.R + + - name: Run R tests + id: rtests + continue-on-error: true + run: Rscript -e 'testthat::test_dir("tests/r")' + + - name: Write result file + if: always() + shell: bash + run: | + mkdir -p result + cat > result/result.json <&2 + exit 1 + fi + echo "version=$VERSION" >> "$GITHUB_OUTPUT" + + - name: Read main version + id: main_version + run: | + git fetch origin main --depth=1 + VERSION=$(git show origin/main:pyproject.toml | grep -m1 '^version *= *' | sed -E 's/version *= *"([^"]+)"/\1/') + if [ -z "$VERSION" ]; then + echo "Could not read version from main's pyproject.toml" >&2 + exit 1 + fi + echo "version=$VERSION" >> "$GITHUB_OUTPUT" + + - name: Install uv + uses: astral-sh/setup-uv@ae62891fec2bb8e7d6c99fc78c9fec3a63790f8d + + - name: Compare versions + run: | + PR_VERSION="${{ steps.pr_version.outputs.version }}" + MAIN_VERSION="${{ steps.main_version.outputs.version }}" + echo "PR version: $PR_VERSION" + echo "main version: $MAIN_VERSION" + uv run --with packaging python3 - "$PR_VERSION" "$MAIN_VERSION" <<'EOF' + import sys + from packaging.version import Version, InvalidVersion + + pr_raw, main_raw = sys.argv[1], sys.argv[2] + + try: + pr = Version(pr_raw) + main = Version(main_raw) + except InvalidVersion as e: + print(f"Could not parse version as PEP 440: {e}", file=sys.stderr) + sys.exit(1) + + if pr <= main: + print(f'Version in pyproject.toml ({pr_raw}) is not greater than main ({main_raw}). Versions must strictly increase.', file=sys.stderr) + sys.exit(1) + + print(f'OK: {pr_raw} > {main_raw}') + EOF \ No newline at end of file diff --git a/.github/workflows/pr_ci_report.py b/.github/workflows/pr_ci_report.py new file mode 100644 index 0000000..2794014 --- /dev/null +++ b/.github/workflows/pr_ci_report.py @@ -0,0 +1,131 @@ +''' +pr_ci_report.py + +Helper script for Pull Request CI Report GitHub Actions workflow +Ingests JSON artifacts containing job results and produces a markdown document with the generated result comment +''' + +# Import external dependencies +import glob, json + +# Define constant icons for pass/partial pass/fail +ICON_PASS = '✅' +ICON_PARTIAL = '⚠️' +ICON_FAIL = '❌' + +# load: given a glob pattern, read all matching JSON files into items and return +def load(pattern): + items = [] + for f in glob.glob(pattern): + with open(f) as fh: + items.append(json.load(fh)) + return items + +# os_key: given an item, return the os_name and os_ver variables (for operating system name/version) +def os_key(item): + return (item.get('os_name', 'Unknown'), item.get('os_ver', 'unknown')) + +# summary_icon: given passed and total jobs, return a summary icon depending on proportion of passed jobs +def summary_icon(passed, total): + if total == 0: + return 'n/a' + if passed == total: + return ICON_PASS + if passed == 0: + return ICON_FAIL + return ICON_PARTIAL + +# load Python and R test results +py = load("results/python/**/result.json") +r = load("results/r/**/result.json") + +# extract all OS keys and sort in descending order +os_keys = sorted(set(os_key(i) for i in py) | set(os_key(i) for i in r), reverse=True) + +# extract all Python versions and sort in ascending order +python_versions = sorted(set(i['python_version'] for i in py), reverse=False) + +# extract all R versions and sort in ascending order +r_versions = sorted(set(i['r_version'] for i in r), reverse=False) + + +# build summary comment as list of lines +lines = [] +lines.append('# Pull Request CI Report') + +# add caution alert about Crux/TRFP tests to comment +lines.append('') +lines.append('> [!CAUTION]') +lines.append('> GitHub Actions runners do not have ThermoRawFileParser or Crux installed. Tests dependent on these binaries are skipped and should instead be run manually before merging.') + +# add summary to comment +lines.append('') +lines.append('
') +lines.append('') +lines.append('## Summary') +lines.append('') +lines.append('Operating System | Source distributions built? | CLI available? | Passed Python tests? | Passed R tests?') +lines.append('--|--|--|--|--') +for name, ver in os_keys: + py_cells = [i for i in py if os_key(i) == (name, ver)] + r_cells = [i for i in r if os_key(i) == (name, ver)] + build_ok = len(py_cells) > 0 and all(c['build'] == 'pass' for c in py_cells) + cli_ok = len(py_cells) > 0 and all(c['cli'] == 'pass' for c in py_cells) + py_passed = sum(1 for c in py_cells if c['tests'] == 'pass') + r_passed = sum(1 for c in r_cells if c['tests'] == 'pass') + build_icon = ICON_PASS if build_ok else ICON_FAIL + cli_icon = ICON_PASS if cli_ok else ICON_FAIL + lines.append(f'{name} ({ver}) | {build_icon} | {cli_icon} | {summary_icon(py_passed, len(py_cells))} | {summary_icon(r_passed, len(r_cells))}') + +# add Python test results to comment +lines.append('') +lines.append('
') +lines.append('') +lines.append('## Python Tests') +lines.append('') +lines.append(f'Python versions tested: {", ".join(python_versions) if python_versions else "none"}') +if python_versions: + lines.append('') + lines.append(' Operating System | ' + ' | '.join(f'Python {v}' for v in python_versions)) + lines.append('-- |' * (len(python_versions) + 1)) + for name, ver in os_keys: + row = [f'{name} ({ver})'] + for v in python_versions: + match = next( + (result for result in py if os_key(result) == (name, ver) and result['python_version'] == v), + None, + ) + row.append('n/a' if match is None else (ICON_PASS if match['tests'] == 'pass' else ICON_FAIL)) + lines.append(' | '.join(row)) + +# add R test results to comment +lines.append('') +lines.append('
') +lines.append('') +lines.append('## R Tests') +lines.append('') +lines.append(f'R versions tested: {", ".join(r_versions) if r_versions else "none"}') +if r_versions: + lines.append('') + lines.append(' Operating System | ' + ' | '.join(f'R {v}' for v in r_versions)) + lines.append('-- |' * (len(r_versions) + 1)) + for name, ver in os_keys: + row = [f'{name} ({ver})'] + for v in r_versions: + match = next( + (result for result in r if os_key(result) == (name, ver) and result['r_version'] == v), + None, + ) + row.append('n/a' if match is None else (ICON_PASS if match['tests'] == 'pass' else ICON_FAIL)) + lines.append(' | '.join(row)) + +# add footer to comment +lines.append('') +lines.append('
') +lines.append('') +lines.append('---') +lines.append('_Auto-generated by the comMS Pull Request CI Report workflow. See full logs in the Actions run for this PR._') + +# write comment to markdown file +with open('pr_ci_report.md', 'w') as f: + f.write('\n'.join(lines)) \ No newline at end of file diff --git a/.github/workflows/publish-release.yml b/.github/workflows/publish-release.yml new file mode 100644 index 0000000..126d199 --- /dev/null +++ b/.github/workflows/publish-release.yml @@ -0,0 +1,40 @@ +name: Create comMS Release + +on: + push: + branches: [main] + +permissions: + contents: write + +jobs: + release: + name: Build and publish release + runs-on: ubuntu-latest + steps: + - uses: actions/checkout@v7 + + - name: Install uv + uses: astral-sh/setup-uv@ae62891fec2bb8e7d6c99fc78c9fec3a63790f8d + + - name: Read version + id: version + run: | + VERSION=$(grep -m1 '^version *= *' pyproject.toml | sed -E 's/version *= *"([^"]+)"/\1/') + if [ -z "$VERSION" ]; then + echo "Could not read version from pyproject.toml" >&2 + exit 1 + fi + echo "version=$VERSION" >> "$GITHUB_OUTPUT" + + - name: Build package + run: uv build + + - name: Publish release + uses: softprops/action-gh-release@v3 + with: + tag_name: v${{ steps.version.outputs.version }} + name: v${{ steps.version.outputs.version }} + generate_release_notes: true + draft: false + files: dist/* \ No newline at end of file diff --git a/CONTRIBUTING.md b/CONTRIBUTING.md index e64f96f..38cb09c 100644 --- a/CONTRIBUTING.md +++ b/CONTRIBUTING.md @@ -31,21 +31,27 @@ See [Installing external tools](./README.md#installing-external-tools) in [`READ comMS follows the below structure: ``` comMS/ - bin/ # External binaries + bin/ # External binaries src/ comms/ - cli/ # Typer command definitions (argument parsing, help text) - commands/ # Logic for each individual command - utils/ # Shared utilities (e.g I/O, paths, config, Crux/TRFP wrappers) - config.toml # Bundled default configuration + cli/ # Typer command definitions (argument parsing, help text) + commands/ # Logic for each individual command + gui/ # PySide6 GUI for the `experiment` command + r/ # R scripts used by the `report` command + deps/ # Dependency checking/installation (also exposed via `comms r-utils`) + sections/ # One script per report section (qc, pca, da, secondary-species, concordance, aux/ev-markers) + utils/ # Shared R helpers (import, normalisation, status tracking, theming) + utils/ # Shared Python utilities (e.g I/O, paths, config, Crux/TRFP wrappers) + config.toml # Bundled default configuration + main.py # Typer app assembly, imported by the `comms` console-script entry point tests/ - conftest.py # Shared fixtures and binary-availability guards + conftest.py # Shared fixtures and binary-availability guards fixtures/ generate_fixtures.py # Synthetic FASTA and mzML generator - unit/ # Pure logic tests, no external binaries required - integration/ # End-to-end tests, may require Crux and/or TRFP - pyproject.toml # Project configuration - uv.lock # UV lockfile + unit/ # Pure logic tests, no external binaries required + integration/ # End-to-end tests, may require Crux and/or TRFP + pyproject.toml # Project configuration + uv.lock # UV lockfile ``` ### CLI vs commands separation diff --git a/README.md b/README.md index 4168a9e..c251cf7 100644 --- a/README.md +++ b/README.md @@ -43,9 +43,14 @@ Tool | Minimum version | Purpose | Platform notes [ThermoRawFileParser][trfp-url] | 1.4.5 | `.RAW` → `.mzML` conversion | Versions < 2.0.0 require [Mono](https://mono-project.com) on Linux/macOS ### `comms report` dependencies -The `report` command requires R (≥ 4.3.0) and a set of R packages (listed in the [report command documentation](./docs/commands.md#the-report-command)). Required R packages can be installed by running: +The `report` command requires R (≥ 4.4.0) and a set of R packages (listed in the [report command documentation](./docs/commands.md#the-report-command)). Check or install the required R packages with: ```bash -Rscript src/comms/r/install_deps.R +comms r-utils check +comms r-utils install +``` +Alternatively, run the underlying script directly: +```bash +Rscript src/comms/r/deps/install_deps.R ``` ## Installation @@ -93,7 +98,7 @@ comms experiment --headless # via terminal # 2. Run comMS analysis pipeline comms pipeline -e /path/to/experiment/dir ``` -Use `--skip-convert` if `.mzML` files are already available, and `--skip-report` to omit the report step. +Use `--skip-convert` if `.mzML` files are already available, and `--skip-report` to omit the report step. `--skip-lfq` and `--skip-quant` omit the two quantification stages individually, `--param-medic` estimates search tolerances before searching, and `--organism-tags`/`-o` supplies per-organism FDR patterns at runtime (see [Per-organism FDR][docs-commands]). Run `comms pipeline --help` for the full option list. You need two inputs to run the pipeline: a sample sheet (TSV or CSV) and a combined FASTA database containing your proteome(s) and contaminants. Both are described in [Input files](./docs/commands.md#input-files). The `--experiment-dir` option sets where comMS reads its configuration and writes its results, explained in [Configuration][docs-config]. diff --git a/docs/commands.md b/docs/commands.md index 0e2e3a7..23d52b0 100644 --- a/docs/commands.md +++ b/docs/commands.md @@ -49,12 +49,12 @@ Command | Description `experiment` | Build a sample sheet, configuration, and metadata (see [Configuration](./configuration.md#creating-an-experiment-the-experiment-command)) `config` | Manage a configuration file (see [Configuration][docs-config]) -### Pipeline +### Pipelines Command | Description ---|--- `pipeline` | Run the full analysis pipeline end-to-end from a sample sheet -### Analysis commands +### Protein Identification Command | Description ---|--- `convert` | Convert `.RAW` files to indexed `.mzML` using ThermoRawFileParser @@ -63,17 +63,23 @@ Command | Description `rescore` | Rescore PSMs using Crux `percolator` on the combined database, then split by organism and run `percolator` per organism for calibrated per-organism FDR `lfq` | Run MS1 label-free quantification using grouped fractions `quantify` | Compute dNSAF spectral counts using Crux `spectral-counts` + +### Downstream Analysis +Command | Description +---|--- `report` | Generate a static analysis report from `quantify` and, optionally, `lfq` output ### Utilities Command | Description ---|--- `license` | Print the comMS license +`r-utils check` | Check whether the R packages required by `report` are installed +`r-utils install` | Install any missing R packages required by `report` `uninstall` | Remove comMS-generated files and print the uninstall command where possible `version` | Print the installed comMS version ## Per-organism FDR -When a combined multi-species FASTA is searched, comMS can apply picked-protein FDR separately per organism. The patterns that split the database can be set once in configuration (`comms config set --organism`, see the [configuration reference](./config-reference.md#protocol-flags)), or supplied at runtime to `rescore` and `pipeline` with `--organism-tags` which takes a comma-separated list of alternating label and pattern pairs: +When a combined multi-species FASTA is searched, comMS can apply picked-protein FDR separately per organism. The patterns that split the database can be set once in configuration (`comms config --organism`, see the [configuration reference](./config-reference.md#protocol-flags)), or supplied at runtime to `rescore` and `pipeline` with `--organism-tags` which takes a comma-separated list of alternating label and pattern pairs: ```bash comms rescore search_dir/ \ --database combined_proteome.fasta \ @@ -87,10 +93,10 @@ comMS then splits the combined FASTA into one sub-FASTA per organism (appending The `report` command runs a set of R-based analysis sections over the quantification output, writing figures and spreadsheets. It requires R (≥ 4.3.0) and the packages listed below. ### Enabling the report command -The report command can be run manually via `comms report` or as part of a `comms pipeline` run. Creating an experiment via `comms experiment` allows selection of reference protein annotation and contaminant files, which are recorded in experiment.toml under [report] and picked up automatically by the `report` command. They can also be supplied at runtime via the `-r`/`--ref-info` and `-c`/`--cont-csv` option flags. +The report command can be run manually via `comms report` or as part of a `comms pipeline` run. Creating an experiment via `comms experiment` allows selection of reference protein annotation and contaminant files, which are recorded in experiment.toml under [report] and picked up automatically by the `report` command. They can also be supplied at runtime via the `-r`/`--ref-info` and `-c`/`--cont-csv` option flags, alongside `-o`/`--organism-prefix`, `-q`/`--quantify-dir`, `-s`/`--sample-sheet`, and `-l`/`--lfq-dir` to override the other resolved inputs. `--min-reps`, `--lfc-threshold`, `--fdr-threshold`, and `--top-n` override the corresponding `[report]` config values for a single run (and are recorded to `report.config.toml` in the output directory when used); `--overwrite` allows writing into an existing report output directory; `--rscript` points to a non-default `Rscript` binary. ### Sections -By default the `qc`, `pca`, and `da` sections run. Use `--section` (repeatable) to choose specific sections, or `--all` to run every section. +By default the `qc`, `pca`, and `da` sections run. Use `--section` (repeatable) to choose specific sections, or `--all` to run every section, including the auxiliary `ev-markers` section. Section | Content ---|--- @@ -99,19 +105,27 @@ Section | Content `da` | limma-based differential abundance per fraction, with volcano plots and DA Venn diagrams `secondary-species` | Secondary-organism proteins per fraction, with a Venn diagram and candidate table `concordance` | LFQ vs dNSAF log₂FC concordance scatter and Venn diagrams (skipped automatically if no `--lfq-dir` is provided) +`ev-markers` | MISEV2023-informed marker-category heatmaps per organism (auxiliary; not run by default) ### Analysis notes **Normalisation.** Normalisation is deliberately not applied across fractions, because the fractions are expected to differ genuinely in protein composition. **Differential abundance.** limma with empirical Bayes shrinkage supports analysis at *n* = 3 ([Ritchie et al., 2015](https://doi.org/10.1093/nar/gkv007), doi:10.1093/nar/gkv007). A Benjamini-Hochberg false-discovery rate is applied within each fraction independently, since across-fraction comparisons are likely confounded by run order. -**EV markers (auxiliary).** A script for analysis of extracellular vesicle marker proteins is provided at `r/sections/aux/ev-markers.R`, informed by the MISEV 2023 guidelines ([Welsh et al., 2024](https://doi.org/10.1002/jev2.12404), doi:10.1002/jev2.12404). It is run manually rather than through `--section`. +**EV markers (auxiliary).** A script for analysis of extracellular vesicle marker proteins is provided at `r/sections/aux/ev-markers.R`, informed by the MISEV 2023 guidelines ([Welsh et al., 2024](https://doi.org/10.1002/jev2.12404), doi:10.1002/jev2.12404). It is not a default section, but `--section ev-markers` can be used to run it on its own, or `--all` to include it alongside every other section. ### R packages -Install the required R packages by running: +Check or install the required R packages via: + +```bash +comms r-utils check +comms r-utils install +``` + +or by running the underlying script directly: ```bash -Rscript src/comms/r/install_deps.R +Rscript src/comms/r/deps/install_deps.R ``` If R is not available, `report` exits with an informative error. @@ -127,6 +141,8 @@ Package | Source `UpSetR` | CRAN `pheatmap` | CRAN `VennDiagram` | CRAN +`iq` | CRAN +`jsonlite` | CRAN `limma` | Bioconductor --- diff --git a/docs/config-reference.md b/docs/config-reference.md index 662791a..7da2d2f 100644 --- a/docs/config-reference.md +++ b/docs/config-reference.md @@ -6,16 +6,20 @@ # Configuration reference -This page lists the values comMS reads from its configuration file: the protocol flags applied by `comms config set`, and the default index, search, and Percolator parameters. For where configuration files live and how they are resolved, see [Configuration][docs-config]. +This page lists the values comMS reads from its configuration file: the protocol flags applied by `comms config`, and the default convert, index, search, rescore, quantify, and report parameters. For where configuration files live and how they are resolved, see [Configuration][docs-config]. ## Contents - [Protocol flags](#protocol-flags) +- [Other configuration flags](#other-configuration-flags) +- [Default convert parameters](#default-convert-parameters) - [Default index parameters](#default-index-parameters) - [Default search parameters](#default-search-parameters) -- [Percolator settings](#percolator-settings) +- [Default rescore parameters](#default-rescore-parameters) +- [Default quantify parameters](#default-quantify-parameters) +- [Default report parameters](#default-report-parameters) ## Protocol flags -`comms config set` applies experiment-specific presets. Flags can be combined in a single call, and most take a positive and a negative form (for example `--ox` / `--no-ox`). +`comms config` applies experiment-specific presets when passed value-setting flags (see [Editing configuration][docs-config] for how the target file is resolved). Flags can be combined in a single call, and most take a positive and a negative form (for example `--ox` / `--no-ox`). This first table covers the peptide-modification and digestion presets used by `index`; the full set of remaining configuration flags, covering `convert`, `search`, `rescore`, `quantify`, and `report`, is in [Other configuration flags](#other-configuration-flags) below. Flag | Effect | Mass | Config key ---|---|---|--- @@ -25,8 +29,9 @@ Flag | Effect | Mass | Config key `--n-cyc` / `--no-n-cyc` | N-terminal Gln to pyro-Glu cyclisation | `1Q-17.027` | `index.nterm_peptide_mods_spec` `--n-ace` / `--no-n-ace` | Protein N-terminal acetylation | `1X+42.011` | `index.nterm_protein_mods_spec` `--clip-met` / `--no-clip-met` | Duplicate peptides with the N-terminal methionine clipped | n/a | `index.clip_n_met` +`--missed-cleavages` | Number of missed enzymatic cleavages allowed | n/a | `index.missed_cleavages` `--high-res` / `--low-res` | Instrument resolution preset | n/a | `search.mz_bin_width`, `search.score_function` -`--custom` | Add or clear a custom Tide `mods_spec` entry | user-defined | `custom_mods` +`--custom` | Add or clear a custom Tide `mods_spec` entry | user-defined | `index.custom_mods` `--organism` | Define label-to-pattern pairs for per-organism FDR | n/a | `[organism]` A few flags need more explanation: @@ -36,18 +41,18 @@ A few flags need more explanation: **Custom modifications (`--custom`)** Custom entries are stored separately and merged with the named-flag modifications at search time. The flag is repeatable, so several entries can be added in one call, and passing an empty string clears them all: ```bash -comms config set --custom "1K+28.0313" # add one entry -comms config set --custom "1K+28.0313" --custom "1R+14.0157" # add several -comms config set --custom "" # clear all custom entries +comms config --custom "1K+28.0313" # add one entry +comms config --custom "1K+28.0313" --custom "1R+14.0157" # add several +comms config --custom "" # clear all custom entries ``` -Passing a modification that is already managed by a named flag (for example `1M+15.9949`, which belongs to `--ox`) produces a warning and is not added. Use the named flag instead. `comms config list` shows both the named-flag and custom values. +Passing a modification that is already managed by a named flag (for example `1M+15.9949`, which belongs to `--ox`) produces a warning and is not added. Use the named flag instead. Running `comms config` with no flags shows both the named-flag and custom values. **Instrument resolution (`--high-res` / `--low-res`)** ```bash -comms config set --high-res # mz_bin_width = 0.02, score_function = xcorr (default) -comms config set --low-res # mz_bin_width = 1.0005079, score_function = combined-p-value +comms config --high-res # mz_bin_width = 0.02, score_function = xcorr (default) +comms config --low-res # mz_bin_width = 1.0005079, score_function = combined-p-value ``` Use `--high-res` for Orbitrap data, the default for modern instruments. Use `--low-res` for ion-trap MS2 data, such as that from older LTQ instruments. @@ -55,11 +60,45 @@ Use `--high-res` for Orbitrap data, the default for modern instruments. Use `--l **Organism patterns (`--organism`)** ```bash -comms config set --organism Org1=Pattern1 Org2=Pattern2 +comms config --organism Org1=Pattern1 Org2=Pattern2 ``` Each argument takes the form `Label=Pattern`, where `Pattern` is matched as a regular expression against FASTA headers. The pairs are used to split a combined FASTA by organism during the rescore step, which enables per-organism picked-protein FDR. Once set, they are applied automatically by `pipeline` and `rescore` unless overridden with `--organism-tags` on the command line. See the [rescore command documentation]((./commands.md#per-organism-fdr)) for the runtime form. +## Other configuration flags +These flags configure `convert`, `search`, `rescore`, `quantify`, and `report` behaviour. As with the protocol flags above, each can be set on a config file via `comms config`, and most can also be overridden for a single run by passing the same flag directly to the corresponding command. + +Flag | Effect | Config key +---|---|--- +`--gzip` / `--no-gzip` | Gzip-compress mzML output | `convert.gzip` +`--format` | ThermoRawFileParser output format code | `convert.format` +`--metadata` | ThermoRawFileParser metadata capture code | `convert.metadata` +`--score-function` | Tide-search score function | `search.score_function` +`--min-peaks` | Minimum peaks required per spectrum | `search.min_peaks` +`--precursor-tolerance-ppm` | Precursor mass tolerance (ppm) | `search.precursor_tolerance_ppm` +`--mz-bin-width` | Fragment m/z bin width (Da) | `search.mz_bin_width` +`--threads` | Default number of threads | `search.threads` +`--protein-enzyme` | Enzyme used for protein-level picked-FDR grouping | `rescore.protein_enzyme` +`--picked-protein` / `--no-picked-protein` | Use picked-protein FDR | `rescore.picked_protein` +`--shared-psm` | Policy for PSMs shared between organisms (`drop`/`include`) | `rescore.shared_psm` +`--measure` | Spectral-counting measure (`NSAF`/`dNSAF`/`SIN`/`EMPAI`) | `quantify.measure` +`--qvalue-threshold` | PSM q-value threshold for quantification | `quantify.qvalue_threshold` +`--unique-mapping` / `--no-unique-mapping` | Require unique peptide-to-protein mapping | `quantify.unique_mapping` +`--min-reps` | Minimum replicates per fraction-treatment group | `report.min_reps` +`--lfc-threshold` | \|log2FC\| threshold for differential abundance | `report.lfc_threshold` +`--fdr-threshold` | BH-FDR threshold for differential abundance | `report.fdr_threshold` +`--top-n` | Number of top DA proteins labelled per volcano plot | `report.top_n_proteins` + +`--score-function`, `--min-peaks`, `--precursor-tolerance-ppm`, and `--mz-bin-width` are also accepted directly by `search`, where `--precursor-tolerance-ppm`/`--mz-bin-width` take priority over `--param-medic` if both are given for the same run. + +## Default convert parameters + +Parameter | Default | Description | Config key +-- | -- | -- | -- +Format | 2 | ThermoRawFileParser output format code (2 = indexed mzML) | `convert.format` +Gzip | True | Compress mzML output | `convert.gzip` +Metadata | 0 | ThermoRawFileParser metadata capture code (0 = JSON metadata) | `convert.metadata` + ## Default index parameters Peptide indices are generated with the following parameters: @@ -83,18 +122,41 @@ Any proteome processed with the `index` command should also include contaminant ## Default search parameters The default configuration applies the following search parameters, informed by [Svozil & Baerenfaller, 2017](https://doi.org/10.1016/bs.mie.2016.11.007) (doi:10.1016/bs.mie.2016.11.007): -Parameter | Default | Description --- | -- | -- -Precursor tolerance | 10 ppm | Precursor mass window -Minimum peaks | 10 | Minimum peaks required per spectrum -Bin width | 0.02 | `mz_bin_width` for high-resolution data (see [instrument resolution](#protocol-flags)) -Score function | xcorr | Scoring for high-resolution data (see [instrument resolution](#protocol-flags)) -Threads | 2 | Default search threads +Parameter | Default | Description | Config key +-- | -- | -- | -- +Precursor tolerance | 10 ppm | Precursor mass window | search.precursor_tolerance_ppm +Minimum peaks | 10 | Minimum peaks required per spectrum | search.min_peaks +Bin width | 0.02 | `mz_bin_width` for high-resolution data (see [instrument resolution](#protocol-flags)) | search.mz_bin_width +Score function | xcorr | Scoring for high-resolution data (see [instrument resolution](#protocol-flags)) | search.score_function +Threads | 2 | Default search threads | search.threads -## Percolator settings +## Default rescore parameters By default, PSM rescoring uses picked-protein FDR ([Savitski et al., 2015](https://doi.org/10.1074/mcp.M114.046995), doi:10.1074/mcp.M114.046995) at a 1% PSM-level FDR threshold, requiring at least two unique peptides per protein for a confident identification. -When a combined multi-species FASTA is used, picked-protein FDR is applied separately per organism. Organism patterns are configured with `comms config set --organism`, or supplied at runtime with `--organism-tags` on the `rescore` and `pipeline` commands. See the [rescore command documentation](./commands.md#per-organism-fdr) for details. For multi-species analyses, the handling policy used in the case of shared PSMs can be set to `'drop'` (default; do not include in either organism) or `'include'` (include in both organisms). The latter option may inflate downstream spectral counts for shared PSMs, so the default option is `'drop'`. +When a combined multi-species FASTA is used, picked-protein FDR is applied separately per organism. Organism patterns are configured with `comms config --organism`, or supplied at runtime with `--organism-tags` on the `rescore` and `pipeline` commands. See the [rescore command documentation](./commands.md#per-organism-fdr) for details. For multi-species analyses, the handling policy used in the case of shared PSMs can be set to `'drop'` (default; do not include in either organism) or `'include'` (include in both organisms). The latter option may inflate downstream spectral counts for shared PSMs, so the default option is `'drop'`. + +Parameter | Default | Description | Config key +-- | -- | -- | -- +Protease | Trypsin | Protease used for protein digestion | rescore.protein_enzyme +Picked-protein FDR | True | Enable picked-protein FDR | rescore.picked_protein +Shared PSM strategy | drop | Strategy used to handle PSMs shared between organisms | rescore.shared_psm + +## Default quantify parameters + +Parameter | Default | Description | Config key +-- | -- | -- | -- +Measure | `dNSAF` | Spectral-counting measure | `quantify.measure` +Q-value threshold | 0.01 | PSM q-value threshold for inclusion | `quantify.qvalue_threshold` +Unique mapping | True | Require unique peptide-to-protein mapping | `quantify.unique_mapping` + +## Default report parameters + +Parameter | Default | Description | Config key +-- | -- | -- | -- +Minimum replicates | 3 | Minimum replicates per fraction-treatment group | `report.min_reps` +Log2FC threshold | 1.0 | \|log2FC\| threshold for differential abundance | `report.lfc_threshold` +FDR threshold | 0.05 | BH-FDR threshold for differential abundance | `report.fdr_threshold` +Top-N proteins | 20 | Number of top DA proteins labelled per volcano plot | `report.top_n_proteins` --- diff --git a/docs/configuration.md b/docs/configuration.md index 64a7432..790c656 100644 --- a/docs/configuration.md +++ b/docs/configuration.md @@ -76,14 +76,14 @@ The simplest way to create an experiment directory is the [`experiment` command] comMS reads its settings from a TOML file. Three sources can supply that file: - **Bundled defaults.** A default `config.toml` ships inside the package and provides the baseline values. It is never edited directly. -- **Global user configuration.** A per-user file created with `comms config init`. It applies to every run unless an experiment provides its own local configuration file. Its location depends on the operating system: +- **Global user configuration.** A per-user file created with `comms config --global`. It applies to every run unless an experiment provides its own local configuration file. Its location depends on the operating system: OS | Path -- | -- Linux/macOS | `~/.config/comms/config.toml` Windows | `%APPDATA%\comms\config.toml` -- **Local configuration.** A `config.toml` inside an experiment's `comms/` folder. It applies only to that experiment. A local file is written by the [`experiment` command](#creating-an-experiment-the-experiment-command), or created with `comms config init -c `. +- **Local configuration.** A `config.toml` inside an experiment's `comms/` folder. It applies only to that experiment. A local file is written by the [`experiment` command](#creating-an-experiment-the-experiment-command), or created by pointing `comms config` at the experiment directory (e.g. `comms config `), which creates the file from defaults if it doesn't already exist. ### Resolution order For any given run, comMS uses the first configuration it finds, in this order: @@ -112,25 +112,40 @@ export COMMS_BIN_DIR=/absolute/path/to/comMS/bin ``` ## Editing configuration: the `config` command -The `config` command edits individual sections of a single configuration file. By default it targets the global user file. Pass `-c` / `--config` with a path to target a local file instead, or `-c global` to be explicit about the global one. +`config` is a single command, not a set of subcommands: it edits or inspects one configuration file, and which behaviour you get depends on which flags you pass. -Subcommand | Purpose +```bash +comms config [PATH] [--global] [--verify | --reset [--force]] [protocol/value flags...] +``` + +Argument/flag | Purpose ---|--- -`init` | Create a configuration file from the defaults (will not overwrite an existing file) -`exists` | Report whether the file exists and print its path -`list` | Print current values, highlighting any that differ from the defaults -`verify` | Check that all expected keys are present -`reset` | Overwrite the file with the defaults (prompts unless `--force`) -`set` | Change values via named flags +`PATH` (positional, optional) | Experiment directory whose local `config.toml` should be targeted. Defaults to the current directory. +`--global` | Target the global user config instead of a local one. +`--verify` | Check that all expected keys are present in the target file (and no unexpected ones), then exit. +`--reset` | Overwrite the target file with comMS defaults (prompts for confirmation unless `--force` is also given), then exit. +`--force` | Skip the confirmation prompt when used with `--reset`. +*(any protocol/value flag, e.g. `--ox`, `--iodo`, `--organism`)* | Apply that value to the target file. + +If none of `--verify`, `--reset`, or a value-setting flag is given (or a value-setting flag is given but resolves to no change), `config` falls back to listing the current values against the defaults — this is also what running `comms config` on its own does. + +**Resolving which file to edit.** Without `PATH` or `--global`, `config` looks in the current directory for a bare `config.toml` or a `config.toml` nested under `comms/`. If neither exists, it offers to create one at `/comms/config.toml`. If a path or `--global` is given and the target file doesn't exist yet, it is created from the bundled defaults first. For example, to view the global configuration, then enable methionine oxidation and cysteine carbamidomethylation in an experiment's local file: ```bash -comms config list -comms config set -c my_experiment/comms/config.toml --ox --iodo +comms config --global +comms config my_experiment --ox --iodo +``` + +To check a config file is complete, or reset it to defaults: + +```bash +comms config my_experiment --verify +comms config my_experiment --reset --force ``` -The full set of `config set` flags, the modifications they apply, and the default parameters are documented in the [configuration reference][docs-config-ref]. +The full set of value-setting flags, the modifications they apply, and the default parameters are documented in the [configuration reference][docs-config-ref]. ## Creating an experiment: the `experiment` command The `experiment` command builds a complete experiment directory, writing the sample sheet, a local `config.toml`, and `experiment.toml` under `/comms/`. It saves manually writing the sample sheet or running `config` for each value. @@ -140,7 +155,7 @@ comms experiment # graphical setup comms experiment --headless # terminal prompts only ``` -Both the GUI and command-line setups walk through naming the experiment and choosing an output directory, defining treatment and fraction groups, importing a directory of `.RAW` / `.mzML` files, assigning each sample to its groups (replicate numbers auto-assign per treatment and fraction and can be overridden), and previewing the sheet before saving. A configuration panel mirrors `comms config set`, so the local `config.toml` it writes uses the same options. The configuration panel also includes a report settings section, where a reference annotation file (TSV/CSV), a contaminant list (CSV), and a primary organism ID prefix can be set. These are written to experiment.toml under [report] and used automatically when comms report is run against the experiment directory. +Both the GUI and command-line setups walk through naming the experiment and choosing an output directory, defining treatment and fraction groups, importing a directory of `.RAW` / `.mzML` files, assigning each sample to its groups (replicate numbers auto-assign per treatment and fraction and can be overridden), and previewing the sheet before saving. A configuration panel mirrors the `comms config` value-setting flags, so the local `config.toml` it writes uses the same options. The configuration panel also includes a report settings section, where a reference annotation file (TSV/CSV), a contaminant list (CSV), and a primary organism ID prefix can be set. These are written to experiment.toml under [report] and used automatically when comms report is run against the experiment directory. ## `config` vs `experiment` The two commands serve different purposes: diff --git a/docs/output-structure.md b/docs/output-structure.md index 032c1be..aad5b86 100644 --- a/docs/output-structure.md +++ b/docs/output-structure.md @@ -32,7 +32,8 @@ All outputs are written under the experiment root, inside `comms/results/`. Each ├─ pca/ # PCA and dendrogram plots ├─ da/ # differential abundance plots and spreadsheet ├─ secondary_species/ # secondary-species plots and spreadsheet - └─ concordance/ # LFQ vs dNSAF concordance (only if --lfq-dir is provided) + ├─ concordance/ # LFQ vs dNSAF concordance (only if --lfq-dir is provided) + └─ ev_markers/ # MISEV2023 marker-category heatmaps (only if run via --section/--all) ``` If a command's output directory already exists, comMS does not overwrite it. Instead it adds an incremental suffix, for example `search-1/`, then `search-2/`, so earlier results are preserved. @@ -50,10 +51,11 @@ Level | Enabled with | Used for `debug` | `-vv` | Detailed program state for debugging: the command invoked, paths scanned, resolved parameters, and per-item detail `progress` | `-v` | Step-wise progress: per-file and per-stage markers `info` | always on | High-level overview: processing counts and final results +`input` | always on (interactive sessions only) | Confirmation prompts and their answers, for example when `config` or `experiment` need to create a file or confirm a destructive action `warn` | always on | A part of a command did not succeed but comMS continued, for example a single item failed, an optional input was missing, or a fallback value was used `error` | always on | A part of a command did not succeed and comMS could not continue -`-vv` implies `-v`. +`-vv` implies `-v`. The `input` level only prompts interactively when running in a terminal (`stdin` is a TTY); otherwise it falls back to the flag's default silently. ### Log files comMS logs to standard output (the terminal you ran the command from) and to a file named after the command (`comms convert`, for example, writes to `convert.log`). The log file is saved alongside that command's other output, so running the same command again does not overwrite an earlier log. diff --git a/pyproject.toml b/pyproject.toml index 5671528..aa5365e 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -1,7 +1,7 @@ [project] name = "comms" -version = "2.0.2" -description = "Add your description here" +version = "3.0.0" +description = "A command line tool for comparative analysis of proteomic data" readme = "README.md" authors = [ { name = "Sam Holland", email = "hol_sam@icloud.com" } @@ -32,9 +32,13 @@ build-backend = "uv_build" [tool.pytest.ini_options] testpaths = ["tests"] -addopts = "--tb=short -ra --cov=comms --cov-report=term-missing" +addopts = "--tb=short -ra --cov=src/comms --cov-report=term-missing" minversion = "7.0" markers = [ "crux: requires the Crux binary to be present under bin/ (auto-skipped if absent)", "trfp: requires ThermoRawFileParser.exe to be present under bin/ (auto-skipped if absent)", ] + +[tool.coverage.run] +source = ["src"] +relative_files = true diff --git a/src/comms/cli/cli.py b/src/comms/cli/cli.py index 52c3cf5..bcb6dbe 100644 --- a/src/comms/cli/cli.py +++ b/src/comms/cli/cli.py @@ -4,7 +4,8 @@ # -- Import external dependencies import logging, typer -from typing import Annotated +from pathlib import Path +from typing import Annotated, Optional # -- Import internal utility functions from comms.utils.settings import initComms @@ -21,6 +22,7 @@ from comms.cli.pipeline import commsPipeline from comms.cli.config import commsConfig from comms.cli.license import commsLicense +from comms.cli.rutils import commsRUtils from comms.cli.uninstall import commsUninstall from comms.cli.version import commsVersion @@ -49,24 +51,29 @@ comms.add_typer(commsLfq) comms.add_typer(commsQuantify) comms.add_typer(commsReport) -comms.add_typer(commsConfig, name='config', help='Manage comMS configuration', rich_help_panel='comMS Configuration') +comms.add_typer(commsConfig) comms.add_typer(commsLicense) +comms.add_typer(commsRUtils, name='r-utils', help='Check or install required R dependencies', rich_help_panel='Utilities') comms.add_typer(commsUninstall) comms.add_typer(commsVersion) # -- Register experiment command @comms.command(rich_help_panel='comMS Configuration') def experiment( + experiment_dir: Annotated[ + Optional[Path], + typer.Argument(help='Existing experiment directory to edit (leave blank to create a new experiment)') + ] = None, headless: Annotated[ bool, typer.Option('--headless', help='Run setup in terminal instead of GUI') ] = False, ): - '''Set up a comMS experiment (sample sheet + config + metadata)''' + '''Set up a comMS experiment (sample sheet + config + metadata), or edit an existing one''' if headless: - experimentFuncs.run_experiment_headless() + experimentFuncs.run_experiment_headless(experiment_dir) else: - experimentFuncs.launch_experiment_gui() + experimentFuncs.launch_experiment_gui(experiment_dir) # ==================== # Top-level callback: --verbose / --debug flags diff --git a/src/comms/cli/config.py b/src/comms/cli/config.py index 9aa941b..2f1630e 100644 --- a/src/comms/cli/config.py +++ b/src/comms/cli/config.py @@ -4,7 +4,8 @@ # -- Import external dependencies import typer -from typing import Annotated, List, Optional +from pathlib import Path +from typing import Annotated, List, Literal, Optional # -- Import internal functions from comms.commands import config as configFuncs @@ -12,100 +13,69 @@ # -- Initialise config Typer class commsConfig = typer.Typer(add_completion=False, invoke_without_command=True) -# Define config file option used in all commands -_CONFIG_OPT = typer.Option( - '-c', '--config', - help="Config file to edit; a path, or 'global' for the user config [default: global]", -) +# -- Define shared path/global parameters +_PATH_ARG = typer.Argument(help='Path to experiment directory [dim](default: ".")[/dim]') +_GLOBAL_OPT = typer.Option('--global', help='Use the global user config instead of a local config.toml') -# -- Define config callback -@commsConfig.callback(invoke_without_command=True) -def config_callback(ctx: typer.Context) -> None: - if ctx.invoked_subcommand is None: - configFuncs.config_exists() - -# -- Define config command: init -@commsConfig.command(rich_help_panel='Config Commands') -def init(config: Annotated[Optional[str], _CONFIG_OPT] = None): - '''Create a user config file with default settings in the OS config directory''' - configFuncs.config_init(config_path=configFuncs._resolveConfigTarget(config)) - -# -- Define config command: exists -@commsConfig.command(rich_help_panel='Config Commands') -def exists(config: Annotated[Optional[str], _CONFIG_OPT] = None): - '''Report whether a user config file exists and print its path''' - configFuncs.config_exists(config_path=configFuncs._resolveConfigTarget(config)) - -# -- Define config command: list -@commsConfig.command(rich_help_panel='Config Commands') -def list(config: Annotated[Optional[str], _CONFIG_OPT] = None): - '''Print current config values, highlighting differences from bundled defaults''' - configFuncs.config_list(config_path=configFuncs._resolveConfigTarget(config)) - -# -- Define config command: verify -@commsConfig.command(rich_help_panel='Config Commands') -def verify(config: Annotated[Optional[str], _CONFIG_OPT] = None): - '''Check that all expected keys are present in the user config file''' - configFuncs.config_verify(config_path=configFuncs._resolveConfigTarget(config)) - -# -- Define config command: reset -@commsConfig.command(rich_help_panel='Config Commands') -def reset( - config: Annotated[Optional[str], _CONFIG_OPT] = None, - force: Annotated[ - bool, - typer.Option('--force', help='Skip confirmation and immediately overwrite config.toml') - ] = False -): - '''Overwrite the user config file with comMS built-in defaults''' - configFuncs.config_reset(config_path=configFuncs._resolveConfigTarget(config), force=force) - -# -- Define config command: set -@commsConfig.command(rich_help_panel='Config Commands') -def set( - config: Annotated[Optional[str], _CONFIG_OPT] = None, - iodo: Annotated[ - Optional[bool], - typer.Option('--iodo/--no-iodo', help='Add (--iodo) or remove (--no-iodo) carbamidomethylation of cysteine as a static modification'), - ] = None, - ox: Annotated[ - Optional[bool], - typer.Option('--ox/--no-ox', help='Add (--ox) or remove (--no-ox) oxidation of methionine as a variable modification'), - ] = None, - phos: Annotated[ - Optional[bool], - typer.Option('--phos/--no-phos', help='Add (--phos) or remove (--no-phos) phosphorylation of serine/threonine/tyrosine as a variable modification'), - ] = None, - n_cyc: Annotated[ - Optional[bool], - typer.Option('--n-cyc/--no-n-cyc', help='Add (--n-cyc) or remove (--no-n-cyc) cyclisation of peptide N-terminal glutamine to pyro-glutamic acid as a variable modification'), - ] = None, - n_ace: Annotated[ - Optional[bool], - typer.Option('--n-ace/--no-n-ace', help='Add (--n-ace) or remove (--no-n-nace) acetylation of protein N-terminal residue as a variable modification'), - ] = None, - custom: Annotated[ - Optional[str], - typer.Option('--custom', help='Add a custom variable modification following Tide mods_spec format; can be passed multiple times; pass empty string "" to remove all custom modifications') - ] = None, - clip_met: Annotated[ - Optional[bool], - typer.Option('--clip-met/--no-clip-met', help="Include (--clip-met) or don't include (--no-clip-met) duplicate N-terminal peptides with clipped N-terminal methionine") - ] = None, - low_res: Annotated[ - Optional[bool], - typer.Option('--low-res/--high-res', help='Set search parameters for low-resolution (--low-res) or high-resolution (--high-res) instruments'), - ] = None, - organism: Annotated[ - Optional[List[str]], - typer.Option('--organism', help='Set organism header patterns for per-organism picked protein FDR [dim](format: OrganismLabel=Pattern)[/dim]'), - ] = None, -): +# -- Define config command (single command, no subcommands) +@commsConfig.command(rich_help_panel='comMS Configuration') +def config( + # -- Global options -- + path: Annotated[Optional[Path], _PATH_ARG] = None, + global_: Annotated[bool, _GLOBAL_OPT] = False, + verify: Annotated[bool, typer.Option('--verify', help='Check that all expected keys are present')] = False, + reset: Annotated[bool, typer.Option('--reset', help='Overwrite with comMS built-in defaults')] = False, + force: Annotated[bool, typer.Option('--force', help='Skip confirmation when using --reset')] = False, + # -- convert command options -- + gzip: Annotated[Optional[bool], typer.Option('--gzip/--no-gzip', help='Compress mzML output using gzip')] = None, + format: Annotated[Optional[int], typer.Option('--format', help='ThermoRawFileParser output format code', min=0, max=4)] = None, + metadata: Annotated[Optional[int], typer.Option('--metadata', help='ThermoRawFileParser metadata capture code', min=0, max=2)] = None, + # -- index command options -- + iodo: Annotated[Optional[bool], typer.Option('--iodo/--no-iodo', help='Add or remove carbamidomethylation of cysteine as a static modification')] = None, + ox: Annotated[Optional[bool], typer.Option('--ox/--no-ox', help='Add or remove oxidation of methionine as a variable modification')] = None, + phos: Annotated[Optional[bool], typer.Option('--phos/--no-phos', help='Add or remove phosphorylation of serine/threonine/tyrosine as a variable modification')] = None, + n_cyc: Annotated[Optional[bool], typer.Option('--n-cyc/--no-n-cyc', help='Add or remove cyclisation of peptide N-terminal glutamine to pyro-glutamic acid')] = None, + n_ace: Annotated[Optional[bool], typer.Option('--n-ace/--no-n-ace', help='Add or remove acetylation of protein N-terminal residue')] = None, + custom: Annotated[Optional[str], typer.Option('--custom', help='Add a custom variable modification (Tide mods_spec format); use "" to clear all custom mods')] = None, + clip_met: Annotated[Optional[bool], typer.Option('--clip-met/--no-clip-met', help='Include or exclude duplicate N-terminal peptides with clipped N-terminal methionine')] = None, + missed_cleavages: Annotated[Optional[int], typer.Option('--missed-cleavages', help='Number of missed enzymatic cleavages allowed', min=0)] = None, + organism: Annotated[Optional[List[str]], typer.Option('--organism', help='Organism header pattern for per-organism picked protein FDR [dim](format: Label=Pattern)[/dim]')] = None, + low_res: Annotated[Optional[bool], typer.Option('--low-res/--high-res', help='Set score_function/mz_bin_width for low- or high-resolution instruments')] = None, + # -- search command options -- + score_function: Annotated[Optional[str], typer.Option('--score-function', help='Tide-search score function')] = None, + min_peaks: Annotated[Optional[int], typer.Option('--min-peaks', help='Minimum peaks required per spectrum', min=1)] = None, + precursor_tolerance_ppm: Annotated[Optional[float], typer.Option('--precursor-tolerance-ppm', help='Precursor mass tolerance in ppm')] = None, + mz_bin_width: Annotated[Optional[float], typer.Option('--mz-bin-width', help='Fragment m/z bin width in Da')] = None, + threads: Annotated[Optional[int], typer.Option('--threads', help='Default number of threads', min=1)] = None, + # -- percolator command options -- + protein_enzyme: Annotated[Optional[str], typer.Option('--protein-enzyme', help='Enzyme used for protein-level picked-FDR grouping')] = None, + picked_protein: Annotated[Optional[bool], typer.Option('--picked-protein/--no-picked-protein', help='Use picked-protein FDR')] = None, + shared_psm: Annotated[Optional[Literal['drop', 'include']], typer.Option('--shared-psm', help='Policy for PSMs shared between organisms')] = None, + # -- quantify command options -- + measure: Annotated[Optional[Literal['NSAF', 'dNSAF', 'SIN', 'EMPAI']], typer.Option('--measure', help='Spectral-counting measure')] = None, + qvalue_threshold: Annotated[Optional[float], typer.Option('--qvalue-threshold', help='PSM q-value threshold for quantification', min=0.0, max=1.0)] = None, + unique_mapping: Annotated[Optional[bool], typer.Option('--unique-mapping/--no-unique-mapping', help='Require unique peptide-to-protein mapping')] = None, + # -- report command options -- + min_reps: Annotated[Optional[int], typer.Option('--min-reps', help='Minimum replicates per fraction-treatment group', min=1)] = None, + lfc_threshold: Annotated[Optional[float], typer.Option('--lfc-threshold', help='|log2FC| threshold for DA', min=0.0)] = None, + fdr_threshold: Annotated[Optional[float], typer.Option('--fdr-threshold', help='BH-FDR threshold for DA', min=0.0, max=1.0)] = None, + top_n: Annotated[Optional[int], typer.Option('--top-n', help='Number of top DA proteins labelled per volcano plot', min=1)] = None, +) -> None: ''' - Set values in user configuration file + View or edit comMS configurations ''' - configFuncs.config_set( - config_path=configFuncs._resolveConfigTarget(config), + if reset: + configFuncs.config_reset(path, global_, force=force) + return + if verify: + configFuncs.config_verify(path, global_) + return + changed = configFuncs.config_set( + path, + global_, + gzip=gzip, + format=format, + metadata=metadata, iodo=iodo, ox=ox, phos=phos, @@ -113,6 +83,24 @@ def set( n_ace=n_ace, custom=custom, clip_met=clip_met, - low_res=low_res, + missed_cleavages=missed_cleavages, organism=organism, - ) \ No newline at end of file + low_res=low_res, + score_function=score_function, + min_peaks=min_peaks, + precursor_tolerance_ppm=precursor_tolerance_ppm, + mz_bin_width=mz_bin_width, + threads=threads, + protein_enzyme=protein_enzyme, + picked_protein=picked_protein, + shared_psm=shared_psm, + measure=measure, + qvalue_threshold=qvalue_threshold, + unique_mapping=unique_mapping, + min_reps=min_reps, + lfc_threshold=lfc_threshold, + fdr_threshold=fdr_threshold, + top_n=top_n, + ) + if not changed: + configFuncs.config_list(path, global_) \ No newline at end of file diff --git a/src/comms/cli/convert.py b/src/comms/cli/convert.py index 1f2df9f..0ccff82 100644 --- a/src/comms/cli/convert.py +++ b/src/comms/cli/convert.py @@ -27,8 +27,16 @@ def convert( ] = Path('.'), gzip: Annotated[ Optional[bool], - typer.Option('--gzip/--no-gzip', help='Gzip-compress mzML output file(s)') + typer.Option('--gzip/--no-gzip', help='Gzip-compress mzML output file(s) [dim][default: config convert.gzip][/dim]') + ] = None, + format: Annotated[ + Optional[int], + typer.Option('--format', help='ThermoRawFileParser output format code [dim][default: config convert.format][/dim]', min=0, max=4) + ] = None, + metadata: Annotated[ + Optional[int], + typer.Option('--metadata', help='ThermoRawFileParser metadata capture code [dim][default: config convert.metadata][/dim]', min=0, max=2) ] = None, ): ctx = ExperimentContext.resolve(experiment_dir) - convertFuncs.run_convert(data, ctx, gzip) \ No newline at end of file + convertFuncs.run_convert(data, ctx, gzip, format, metadata) \ No newline at end of file diff --git a/src/comms/cli/index.py b/src/comms/cli/index.py index 97c2700..c568a40 100644 --- a/src/comms/cli/index.py +++ b/src/comms/cli/index.py @@ -25,6 +25,25 @@ def index( Optional[Path], typer.Option('-e', '--experiment-dir', help='Experiment directory', exists=True, file_okay=False, dir_okay=True, writable=True) ] = Path('.'), + iodo: Annotated[Optional[bool], typer.Option('--iodo/--no-iodo', help='Override carbamidomethylation of cysteine for this run only [dim][default: config][/dim]')] = None, + ox: Annotated[Optional[bool], typer.Option('--ox/--no-ox', help='Override oxidation of methionine for this run only [dim][default: config][/dim]')] = None, + phos: Annotated[Optional[bool], typer.Option('--phos/--no-phos', help='Override phosphorylation of S/T/Y for this run only [dim][default: config][/dim]')] = None, + n_cyc: Annotated[Optional[bool], typer.Option('--n-cyc/--no-n-cyc', help='Override N-terminal Gln cyclisation for this run only [dim][default: config][/dim]')] = None, + n_ace: Annotated[Optional[bool], typer.Option('--n-ace/--no-n-ace', help='Override N-terminal protein acetylation for this run only [dim][default: config][/dim]')] = None, + custom: Annotated[Optional[str], typer.Option('--custom', help='Add a custom variable modification for this run only [dim][default: config][/dim]')] = None, + clip_met: Annotated[Optional[bool], typer.Option('--clip-met/--no-clip-met', help='Override clipped N-terminal methionine handling for this run only [dim][default: config][/dim]')] = None, + missed_cleavages: Annotated[Optional[int], typer.Option('--missed-cleavages', help='Missed cleavages for this run only [dim][default: config index.missed_cleavages][/dim]', min=0)] = None, ): ctx = ExperimentContext.resolve(experiment_dir) - indexFuncs.run_index(database, ctx) \ No newline at end of file + indexFuncs.run_index( + database, + ctx, + iodo=iodo, + ox=ox, + phos=phos, + n_cyc=n_cyc, + n_ace=n_ace, + custom=custom, + clip_met=clip_met, + missed_cleavages=missed_cleavages, + ) \ No newline at end of file diff --git a/src/comms/cli/quantify.py b/src/comms/cli/quantify.py index 343c1ee..1a8bf76 100644 --- a/src/comms/cli/quantify.py +++ b/src/comms/cli/quantify.py @@ -5,7 +5,7 @@ # -- Import external dependencies import typer from pathlib import Path -from typing import Annotated, Optional +from typing import Annotated, Literal, Optional # -- Import internal functions from comms.commands import quantify as quantifyFuncs @@ -29,6 +29,25 @@ def quantify( Optional[Path], typer.Option('-e', '--experiment-dir', help='Experiment directory', exists=True, file_okay=False, dir_okay=True, writable=True) ] = Path('.'), + measure: Annotated[ + Optional[Literal['NSAF', 'dNSAF', 'EMPAI', 'SIN']], + typer.Option('--measure', help='Spectral-counting measure [dim][default: config quantify.measure][/dim]') + ] = None, + qvalue_threshold: Annotated[ + Optional[float], + typer.Option('--qvalue-threshold', help='PSM q-value threshold for inclusion [dim][default: config quantify.qvalue_threshold][/dim]', min=0.0, max=1.0) + ] = None, + unique_mapping: Annotated[ + Optional[bool], + typer.Option('--unique-mapping/--no-unique-mapping', help='Require unique peptide-to-protein mapping [dim][default: config quantify.unique_mapping][/dim]') + ] = None, ): ctx = ExperimentContext.resolve(experiment_dir) - quantifyFuncs.run_quantify(psm_dir, database, ctx) \ No newline at end of file + quantifyFuncs.run_quantify( + psm_dir, + database, + ctx, + measure=measure, + qvalue_threshold=qvalue_threshold, + unique_mapping=unique_mapping, + ) \ No newline at end of file diff --git a/src/comms/cli/report.py b/src/comms/cli/report.py index a48f738..04a99b8 100644 --- a/src/comms/cli/report.py +++ b/src/comms/cli/report.py @@ -23,15 +23,15 @@ def report( organism_prefix: Annotated[ Optional[str], - typer.Option('-o', '--organism-prefix', help='ID prefix for the primary organism [dim][default: experiment organism prefix][/dim]') + typer.Option('-o', '--organism-prefix', help='ID prefix for the primary organism [dim]\\[default: experiment organism prefix][/dim]') ] = None, quantify_dir: Annotated[ Optional[Path], - typer.Option('-q', '--quantify-dir', help='Path to quantification results [dim][default: quantify output][/dim]') + typer.Option('-q', '--quantify-dir', help='Path to quantification results [dim]\\[default: quantify output][/dim]') ] = None, sample_sheet: Annotated[ Optional[Path], - typer.Option('-s', '--sample-sheet', help='Path to sample sheet [dim][default: experiment sample sheet][/dim]') + typer.Option('-s', '--sample-sheet', help='Path to sample sheet [dim]\\[default: experiment sample sheet][/dim]') ] = None, experiment_dir: Annotated[ Optional[Path], @@ -39,28 +39,32 @@ def report( ] = Path('.'), lfq_dir: Annotated[ Optional[Path], - typer.Option('-l', '--lfq-dir', help='Path to LFQ results [dim][default: lfq output][/dim]') + typer.Option('-l', '--lfq-dir', help='Path to LFQ results [dim]\\[default: lfq output][/dim]') ] = None, ref_info: Annotated[ Optional[Path], - typer.Option('-r', '--ref-info', help='Protein metadata TSV [dim][default: experiment ref_info][/dim]') + typer.Option('-r', '--ref-info', help='Protein metadata TSV [dim]\\[default: experiment ref_info][/dim]') ] = None, cont_csv: Annotated[ Optional[Path], - typer.Option('-c', '--cont-csv', help='Contaminant annotations CSV [dim][default: experiment cont_csv][/dim]') + typer.Option('-c', '--cont-csv', help='Contaminant annotations CSV [dim]\\[default: experiment cont_csv][/dim]') ] = None, min_reps: Annotated[ - int, - typer.Option('--min-reps', help='Minimum replicates per fraction-treatment group', min=1) - ] = 3, + Optional[int], + typer.Option('--min-reps', help='Minimum replicates per fraction-treatment group [dim]\\[default: config report.min_reps][/dim]', min=1) + ] = None, lfc_threshold: Annotated[ - float, - typer.Option('--lfc-threshold', help='|log2FC| threshold for DA', min=0.0) - ] = 1.0, + Optional[float], + typer.Option('--lfc-threshold', help='|log2FC| threshold for DA [dim]\\[default: config report.lfc_threshold][/dim]', min=0.0) + ] = None, fdr_threshold: Annotated[ - float, - typer.Option('--fdr-threshold', help='BH-FDR threshold for DA', min=0.0, max=1.0) - ] = 0.05, + Optional[float], + typer.Option('--fdr-threshold', help='BH-FDR threshold for DA [dim]\\[default: config report.fdr_threshold][/dim]', min=0.0, max=1.0) + ] = None, + top_n: Annotated[ + Optional[int], + typer.Option('--top-n', help='Number of top DA proteins labelled per volcano plot [dim]\\[default: config report.top_n_proteins][/dim]', min=1) + ] = None, section: Annotated[ Optional[list[str]], typer.Option('--section', help='Section(s) to run (repeatable)') @@ -91,6 +95,7 @@ def report( min_reps=min_reps, lfc_threshold=lfc_threshold, fdr_threshold=fdr_threshold, + top_n=top_n, sections=sections, overwrite=overwrite, rscript=rscript, diff --git a/src/comms/cli/rescore.py b/src/comms/cli/rescore.py index d4793ca..f0d273c 100644 --- a/src/comms/cli/rescore.py +++ b/src/comms/cli/rescore.py @@ -5,7 +5,7 @@ # -- Import external dependencies import typer from pathlib import Path -from typing import Annotated, Optional +from typing import Annotated, Literal, Optional # -- Import internal functions from comms.commands import rescore as rescoreFuncs @@ -33,6 +33,26 @@ def rescore( Optional[Path], typer.Option('-e', '--experiment-dir', help='Experiment directory', exists=True, file_okay=False, dir_okay=True, writable=True) ] = Path('.'), + protein_enzyme: Annotated[ + Optional[str], + typer.Option('--protein-enzyme', help='Enzyme used for protein-level picked-FDR grouping [dim][default: config rescore.protein_enzyme][/dim]') + ] = None, + picked_protein: Annotated[ + Optional[bool], + typer.Option('--picked-protein/--no-picked-protein', help='Use picked-protein FDR [dim][default: config rescore.picked_protein][/dim]') + ] = None, + shared_psm: Annotated[ + Optional[Literal['drop', 'include']], + typer.Option('--shared-psm', help='Policy for PSMs shared between organisms [dim][default: config rescore.shared_psm][/dim]') + ] = None, ): ctx = ExperimentContext.resolve(experiment_dir) - rescoreFuncs.run_rescore(psm_dir, database, ctx, organism_tags) \ No newline at end of file + rescoreFuncs.run_rescore( + psm_dir, + database, + ctx, + organism_tags, + protein_enzyme=protein_enzyme, + picked_protein=picked_protein, + shared_psm=shared_psm, + ) \ No newline at end of file diff --git a/src/comms/cli/rutils.py b/src/comms/cli/rutils.py new file mode 100644 index 0000000..616590c --- /dev/null +++ b/src/comms/cli/rutils.py @@ -0,0 +1,31 @@ +''' +comMS CLI subcommand for managing R dependencies +''' + +# -- Import external dependencies +import typer +from typing import Annotated, List, Optional + +# -- Import internal functions +from comms.utils import installrdeps as rUtils +from comms.utils.log import logMsg + +# -- Initialise Typer class +commsRUtils = typer.Typer(add_completion=False, invoke_without_command=True) + +# -- Define rUtils callback +@commsRUtils.callback(invoke_without_command=False) +def rutils_callback(ctx: typer.Context) -> None: + logMsg('r-utils') + +# -- Define R utility command: check +@commsRUtils.command(rich_help_panel='R Utilities') +def check(): + '''Check that all required R dependencies are installed''' + rUtils.check_r_dependencies() + +# -- Define R utility command: install +@commsRUtils.command(rich_help_panel='R Utilities') +def install(): + '''Install any missing R dependencies''' + rUtils.install_r_dependencies_terminal() \ No newline at end of file diff --git a/src/comms/cli/search.py b/src/comms/cli/search.py index d30a5c7..585e67f 100644 --- a/src/comms/cli/search.py +++ b/src/comms/cli/search.py @@ -37,6 +37,32 @@ def search( int, typer.Option('--threads', help='Number of threads', min=1) ] = None, + score_function: Annotated[ + Optional[str], + typer.Option('--score-function', help='Tide-search score function [dim][default: config search.score_function][/dim]') + ] = None, + min_peaks: Annotated[ + Optional[int], + typer.Option('--min-peaks', help='Minimum peaks required per spectrum [dim][default: config search.min_peaks][/dim]', min=1) + ] = None, + precursor_tolerance_ppm: Annotated[ + Optional[float], + typer.Option('--precursor-tolerance-ppm', help='Precursor mass tolerance in ppm; takes priority over --param-medic if both given [dim][default: config search.precursor_tolerance_ppm][/dim]') + ] = None, + mz_bin_width: Annotated[ + Optional[float], + typer.Option('--mz-bin-width', help='Fragment m/z bin width in Da; takes priority over --param-medic if both given [dim][default: config search.mz_bin_width][/dim]') + ] = None, ): ctx = ExperimentContext.resolve(experiment_dir) - searchFuncs.run_search(data, index, ctx, param_medic, threads) \ No newline at end of file + searchFuncs.run_search( + data, + index, + ctx, + param_medic, + threads, + score_function=score_function, + min_peaks=min_peaks, + precursor_tolerance_ppm=precursor_tolerance_ppm, + mz_bin_width=mz_bin_width, + ) \ No newline at end of file diff --git a/src/comms/commands/config.py b/src/comms/commands/config.py index 4768185..2df29e1 100644 --- a/src/comms/commands/config.py +++ b/src/comms/commands/config.py @@ -11,112 +11,139 @@ from rich import print from rich.console import Console from rich.table import Table -from typing import Annotated # -- Import internal functions +from comms.utils.context import _normalise_dirs from comms.utils.log import logMsg -from comms.utils.settings import loadDefaultConfig, globalConfigPath +from comms.utils.modspec import apply_custom_mod, apply_organism, apply_protocol_flags, parse_organism_arg +from comms.utils.settings import loadDefaultConfig, globalConfigPath, _writeConfigTo -# -- Define modification constants -CARBAMIDOMETHYL_MOD = 'C+57.0215' # static carbamidomethylation of Cys -MET_OX_MOD = '1M+15.9949' # variable Met oxidation -PHOSPHO_MOD = '1STY+79.966331' # variable STY phosphorylation -NCYC_MOD = '1Q-17.027' # N-terminal Gln cyclisation -NACE_MOD = '1X+42.011' # N-terminal protein acetylation -MANAGED_MOD_PATTERNS: dict[str, str] = { - r'^\d*C[+\-]': '--iodo / --no-iodo', - r'^\d*M\+15\.9949': '--ox / --no-ox', - r'^\d*STY\+79\.966331': '--phos / --no-phos', -} # mods that --custom is not allowed to duplicate (maps the residue/pattern that identifies each managed mod to its flag name) +# -- _confirm: yes/no prompt via logMsg.input, returned as a bool +def _confirm(msg: str, default: bool) -> bool: + msg = f'{msg} [dim]({"Y/n" if default else "y/N"})[/dim]' + answer = logMsg.input(msg, choices=['y', 'n'], default='y' if default else 'n', case_sensitive=False, show_choices=False, show_default=False) + return str(answer).strip().lower() == 'y' -# -- Define resolution constants -MZ_BIN_WIDTH_HIGH_RES = 0.02 # high-resolution instruments (default) -MZ_BIN_WIDTH_LOW_RES = 1.0005079 # low-resolution instruments -SCORE_FUNC_HIGH_RES = 'xcorr' # high-resolution instruments (default) -SCORE_FUNC_LOW_RES = 'combined-p-value' # low-resolution instruments +# -- _loadConfigFile: returns the config as a dict +def _loadConfigFile(config_path: Path) -> dict: + with config_path.open('rb') as f: + return tomllib.load(f) +# -- _flatten: returns a flat dict from a nested dict, with dot-separated keys +def _flatten(d: dict, prefix: str = '') -> dict: + out = {} + for k, v in d.items(): + key = f'{prefix}.{k}' if prefix else k + if isinstance(v, dict): + out.update(_flatten(v, key)) + else: + out[key] = v + return out -# ========================= # -# DEFINE CONFIG SUBCOMMANDS # -# ========================= # -# -- config_init: creates a config file with default settings in the OS config directory -def config_init(config_path: Path | None = None): - logMsg('config') - config_path = config_path or globalConfigPath() - logMsg.debug(f'Checking config path: {config_path}') - if not _configCheck(config_path, exists=False): - raise SystemExit(1) - try: - logMsg.progress(f'Writing default config to {config_path}') - _writeConfigTo(loadDefaultConfig(), config_path) - logMsg.info(f'Config file written to {config_path}') - except Exception as e: - logMsg.error(f'Failed to write config: {e}') - raise SystemExit(1) +# -- _printTable: prints a Rich table comparing current and default config values +def _printTable(user_config: dict, default_config: dict) -> None: + console = Console() + table = Table(title='comMS configuration', show_header=True, header_style='bold', show_lines=False) + table.add_column('Key', style='cyan', no_wrap=True) + table.add_column('Current value', justify='right') + table.add_column('Default value', justify='right', style='dim') + table.add_column('', width=2) + for key in sorted(default_config.keys()): + default_val = default_config[key] + user_val = user_config.get(key, '[bold red]MISSING[/bold red]') + changed = str(user_val) != str(default_val) + status = '[yellow]≠[/yellow]' if changed else '[green]✓[/green]' + user_str = f'[yellow]{user_val}[/yellow]' if changed else str(user_val) + table.add_row(key, user_str, str(default_val), status) + console.print(table) -# -- config_exists: reports whether a config file exists and prints its path -def config_exists(config_path: Path | None = None): - logMsg('config') - config_path = config_path or globalConfigPath() - logMsg.debug(f'Checking for config at {config_path}') - if config_path.exists(): - logMsg.info(f'Config file found at {config_path}') +# -- _print_diff_summary: prints only the keys that changed between two flattened config dicts +def _print_diff_summary(before: dict, after: dict) -> None: + changed = {k: (before.get(k), v) for k, v in after.items() if before.get(k) != v} + if not changed: + print('\n[dim]No changes made.[/dim]\n') + return + print() + for key, (old, new) in sorted(changed.items()): + print(f'[bold green]✓[/bold green] [dim]{key}[/dim]: [dim]{old}[/dim] → [cyan]{new}[/cyan]') + print() + +# -- _resolve_or_create: returns the Path to edit, creating it from defaults first if needed +def _resolve_or_create(path: Path | None, use_global: bool) -> Path: + ''' + Resolve the config.toml target and make sure it exists, creating it from bundled defaults if not + ''' + if use_global: + target = globalConfigPath() + elif path is not None: + _, comms_dir = _normalise_dirs(path) + target = comms_dir / 'config.toml' else: - logMsg.error(f'No config file at {config_path}') - raise SystemExit(1) + bare, nested = Path.cwd() / 'config.toml', Path.cwd() / 'comms' / 'config.toml' + if bare.exists() and nested.exists(): + logMsg.error(f'Both {bare} and {nested} exist in the current directory. Remove one before running comms config here.') + raise SystemExit(1) + if bare.exists(): + target = bare + elif nested.exists(): + target = nested + else: + logMsg.warn(f'No local config found in the current directory. Did you mean to use [bold]--global[/bold]?') + create_answer = _confirm(msg=f'Create default config at {nested}', default=True) + if not create_answer: + raise SystemExit(0) + target = nested + if not target.exists(): + logMsg.debug(f'Creating default config at {target}') + _writeConfigTo(loadDefaultConfig(), target) + return target # -- config_list: prints current config values, highlighting differences from bundled defaults -def config_list(config_path: Path | None = None): +def config_list(path, global_) -> None: logMsg('config') + logMsg.debug(f'Resolving configuration file') + config_path = _resolve_or_create(path, global_) logMsg.debug(f'Listing config values') - config_path = config_path or globalConfigPath() default_config = _flatten(loadDefaultConfig()) - if _configCheck(config_path, exists=True): - print(f'[bold blue]Current config:[/bold blue] [cyan]{config_path}[/cyan]\n') - current_config = _flatten(_loadConfigFile(config_path)) - else: - print(f'[bold blue]Current config:[/bold blue] built-in defaults\n') - current_config = default_config + print(f'\n[bold blue]Current config:[/bold blue] [cyan]{config_path}[/cyan]\n') + current_config = _flatten(_loadConfigFile(config_path)) _printTable(current_config, default_config) print() # -- config_verify: checks that all expected keys are present in the config file -def config_verify(config_path: Path | None = None): +def config_verify(path, global_) -> None: logMsg('config') - config_path = config_path or globalConfigPath() + logMsg.debug(f'Resolving configuration file') + config_path = _resolve_or_create(path, global_) logMsg.debug(f'Verifying config keys at {config_path}') - if not _configCheck(config_path, exists=True): - logMsg.error(f'No config to verify at {config_path}') - raise SystemExit(1) user_config = _flatten(_loadConfigFile(config_path)) default_config = _flatten(loadDefaultConfig()) missing = [k for k in default_config if k not in user_config] unexpected = [k for k in user_config if k not in default_config] if not missing and not unexpected: - logMsg.info(f'User config {config_path} is valid') + logMsg.info(f'Config {config_path} is valid') return - logMsg.error(f'User config invalid: {len(missing)} missing, {len(unexpected)} unexpected key(s)') + logMsg.error(f'Config invalid: {len(missing)} missing, {len(unexpected)} unexpected key(s)') if missing: - logMsg.warn(f'Missing keys in config: {missing}') print(f'[bold red]ERROR:[/bold red] {len(missing)} missing key(s):') for k in sorted(missing): print(f'\t[red]✗[/red] {k} [dim](expected: {default_config[k]})[/dim]') if unexpected: - logMsg.warn(f'Unexpected keys in config: {unexpected}') print(f'[bold red]ERROR:[/bold red] {len(unexpected)} unexpected key(s):') for k in sorted(unexpected): print(f'\t[red]?[/red] {k}: {user_config[k]}') - print(f'Run [bold]comms config reset[/bold] to restore defaults.\n') + print(f'Run [bold]comms config --reset[/bold] to restore defaults.\n') raise SystemExit(1) # -- config_reset: overwrites the config file with comMS built-in defaults -def config_reset(config_path: Path | None = None, force: bool = False): +def config_reset(path, global_, force: bool = False) -> None: logMsg('config') - config_path = config_path or globalConfigPath() + logMsg.debug(f'Resolving configuration file') + config_path = _resolve_or_create(path, global_) if not force: logMsg.warn(f'This will overwrite {config_path} with comMS defaults.') - if not typer.confirm('All custom settings will be lost. Continue?'): - logMsg.debug(f'Reset cancelled') + if not _confirm('Continue with reset'): + logMsg.debug('Reset cancelled') raise SystemExit(0) try: _writeConfigTo(loadDefaultConfig(), config_path) @@ -125,339 +152,63 @@ def config_reset(config_path: Path | None = None, force: bool = False): logMsg.error(f'Failed to reset config: {e}') raise SystemExit(1) -# -- config_set: apply named flags to the config -def config_set( - config_path: Path | None = None, - iodo: bool | None = None, - low_res: bool | None = None, - organism: list[str] | None = None, - ox: bool | None = None, - phos: bool | None = None, - n_cyc: bool | None = None, - n_ace: bool | None = None, - custom: str | None = None, - clip_met: bool | None = None, -) -> None: - # Set up logger +# -- config_set: apply any given flags to the config file; returns True if anything changed +def config_set(path, global_, **flags) -> bool: logMsg('config') - logMsg.debug(f'Applying set flags: iodo={iodo}; ox={ox}; phos={phos}; n_cyc={n_cyc}; n_ace={n_ace}; low_res={low_res}; organism={organism}; custom={custom!r}; clip_met={clip_met}') - # Check at least one flag set - if all(v is None for v in (iodo, ox, phos, n_cyc, n_ace, low_res, organism, custom, clip_met)): - logMsg.error(f'No flags supplied to config set') - raise SystemExit(1) - # Check if config exists - config_path = config_path or globalConfigPath() - if not config_path.exists(): - logMsg.debug(f'No config found, creating from defaults at {config_path}') - _writeConfigTo(loadDefaultConfig(), path=config_path) - # Load config + logMsg.debug(f'Resolving configuration file') + config_path = _resolve_or_create(path, global_) + logMsg.debug(f'Applying flags: {flags}') + if all(v is None for v in flags.values()): + return False try: cfg = _loadConfigFile(config_path) except Exception as e: logMsg.error(f'Failed to read config file: {e}') raise SystemExit(1) - # Apply any passed flags - cfg = _apply_protocol_flags( + before = _flatten(cfg).copy() + cfg = apply_protocol_flags( cfg, - iodo=iodo, - ox=ox, - phos=phos, - n_cyc=n_cyc, - n_ace=n_ace, - low_res=low_res, - clip_met=clip_met + iodo=flags.get('iodo'), + ox=flags.get('ox'), + phos=flags.get('phos'), + n_cyc=flags.get('n_cyc'), + n_ace=flags.get('n_ace'), + clip_met=flags.get('clip_met'), + low_res=flags.get('low_res'), + missed_cleavages=flags.get('missed_cleavages'), ) - if organism is not None: - cfg = _apply_organism(cfg, _parse_organism_arg(organism)) - if custom is not None: + if flags.get('organism') is not None: + cfg = apply_organism(cfg, parse_organism_arg(flags['organism'])) + if flags.get('custom') is not None: current = cfg.get('index', {}).get('custom_mods', '') - cfg.setdefault('index', {})['custom_mods'] = _apply_custom_mod(current, custom) - # Write updated config + cfg.setdefault('index', {})['custom_mods'] = apply_custom_mod(current, flags['custom']) + direct = { + ('convert', 'gzip'): flags.get('gzip'), + ('convert', 'format'): flags.get('format'), + ('convert', 'metadata'): flags.get('metadata'), + ('search', 'score_function'): flags.get('score_function'), + ('search', 'min_peaks'): flags.get('min_peaks'), + ('search', 'precursor_tolerance_ppm'): flags.get('precursor_tolerance_ppm'), + ('search', 'mz_bin_width'): flags.get('mz_bin_width'), + ('search', 'threads'): flags.get('threads'), + ('rescore', 'protein_enzyme'): flags.get('protein_enzyme'), + ('rescore', 'picked_protein'): flags.get('picked_protein'), + ('rescore', 'shared_psm'): flags.get('shared_psm'), + ('quantify', 'measure'): flags.get('measure'), + ('quantify', 'qvalue_threshold'): flags.get('qvalue_threshold'), + ('quantify', 'unique_mapping'): flags.get('unique_mapping'), + ('report', 'min_reps'): flags.get('min_reps'), + ('report', 'lfc_threshold'): flags.get('lfc_threshold'), + ('report', 'fdr_threshold'): flags.get('fdr_threshold'), + ('report', 'top_n_proteins'): flags.get('top_n'), + } + for (section, key), value in direct.items(): + if value is not None: + cfg.setdefault(section, {})[key] = value try: _writeConfigTo(cfg, config_path) except Exception as e: logMsg.error(f'Failed to write config file: {e}') raise SystemExit(1) - # Print summary - _printSetSummary(iodo=iodo, ox=ox, phos=phos, n_cyc=n_cyc, n_ace=n_ace, low_res=low_res, organism=organism, custom=custom, clip_met=clip_met) - print() - - -# ======================= # -# DEFINE INTERNAL HELPERS # -# ======================= # -# -- _resolveConfigTarget: returns the Path to edit (global user config or a local file) -def _resolveConfigTarget(target: str | None) -> Path: - if target is None or target.upper() == 'GLOBAL': - return globalConfigPath() - return Path(target) - -# -- _loadConfigFile: returns the config as a dict -def _loadConfigFile(config_path: Path | None = None) -> dict: - config_path = config_path or globalConfigPath() - if not config_path.exists(): - raise FileNotFoundError(f'No config found at {config_path}.') - with config_path.open('rb') as f: - return tomllib.load(f) - -# -- _writeConfigTo: writes a config dict to a given path -def _writeConfigTo(config: dict, path: Path) -> None: - path.parent.mkdir(parents=True, exist_ok=True) - with path.open('wb') as f: - tomli_w.dump(config, f) - -# -- _writeConfig: writes config dict to the global config path -def _writeConfig(config: dict): - _writeConfigTo(config, globalConfigPath()) - -# -- _flatten: returns a flat dict from a nested dict, with dot-separated keys -def _flatten(d: dict, prefix: str = '') -> dict: - out = {} - for k, v in d.items(): - key = f'{prefix}.{k}' if prefix else k - if isinstance(v, dict): - out.update(_flatten(v, key)) - else: - out[key] = v - return out - -# -- _configCheck: returns True if the config file existence matches the expected state -def _configCheck(config_path: Path, exists: bool) -> bool: - if exists: - if config_path.exists(): - return True - print(f'\n[bold yellow]WARNING:[/bold yellow] No config at [cyan]{config_path}[/cyan]\nRun [bold]comms config init[/bold] to create one.\n') - return False - else: - if config_path.exists(): - print(f'\n[bold yellow]WARNING:[/bold yellow] Config already exists at [cyan]{config_path}[/cyan]\nRun [bold]comms config reset[/bold] to reset to defaults.\n') - return False - return True - -# -- _printTable: prints a Rich table comparing current and default config values -def _printTable(user_config: dict, default_config: dict): - console = Console() - table = Table(title='comMS configuration', show_header=True, header_style='bold', show_lines=False) - table.add_column('Key', style='cyan', no_wrap=True) - table.add_column('Current value', justify='right') - table.add_column('Default value', justify='right', style='dim') - table.add_column('', width=2) - for key in sorted(default_config.keys()): - default_val = default_config[key] - user_val = user_config.get(key, '[bold red]MISSING[/bold red]') - changed = str(user_val) != str(default_val) - status = '[yellow]≠[/yellow]' if changed else '[green]✓[/green]' - user_str = f'[yellow]{user_val}[/yellow]' if changed else str(user_val) - table.add_row(key, user_str, str(default_val), status) - console.print(table) - - -# ========================= # -# DEFINE CONFIG SET HELPERS # -# ========================= # -# -- _apply_protocol_flags: returns dictionary of config options -def _apply_protocol_flags( - cfg: dict, - *, - iodo: bool | None = None, - ox: bool | None = None, - phos: bool | None = None, - n_cyc: bool | None = None, - n_ace: bool | None = None, - clip_met: bool | None = None, - low_res: bool | None = None, -) -> dict: - ''' - Apply protocol flags to a config dictionary and return it - iodo — owns the Cys slot in index.fixed_mods exclusivel - ox — adds/removes 1M+15.9949 in index.mods_spec - phos — adds/removes 1STY+79.966331 in index.mods_spec - n_cyc — adds/removes 1Q-17.027 in index.nterm_peptide_mods_spec - n_ace — adds/removes 1X+42.011 in index.nterm_protein_mod_spec - low_res — sets search.mz_bin_width and index.score_function - ''' - cfg.setdefault('search', {}) - cfg['index'].setdefault('fixed_mods', '') - cfg['index'].setdefault('nterm_peptide_mods_spec', '') - cfg['index'].setdefault('nterm_protein_mods_spec', '') - if iodo is not None: - cfg['index']['fixed_mods'] = _apply_iodo(cfg['index'].get('fixed_mods', ''), iodo=iodo) - logMsg.debug(f'{'--iodo' if iodo else '--no-iodo'} applied: fixed_mods updated to {cfg['index']['fixed_mods']}') - if ox is not None: - spec = cfg['index'].get('mods_spec', '') - if ox: - cfg['index']['mods_spec'] = _apply_mod(spec, mod=MET_OX_MOD) - else: - cfg['index']['mods_spec'] = _apply_mod(spec, mod='', exclusive_pattern=r'^\d*M\+15\.9949') - logMsg.debug(f'{'--ox' if ox else '--no-ox'} applied: mods_spec updated to {cfg['index']['mods_spec']}') - if phos is not None: - spec = cfg['index'].get('mods_spec', '') - if phos: - cfg['index']['mods_spec'] = _apply_mod(spec, mod=PHOSPHO_MOD) - else: - cfg['index']['mods_spec'] = _apply_mod(spec, mod='', exclusive_pattern=r'^\d*STY\+79\.966331') - logMsg.debug(f'{'--phos' if phos else '--no-phos'} applied: mods_spec updated to {cfg['index']['mods_spec']}') - if n_cyc is not None: - spec = cfg['index'].get('nterm_peptide_mods_spec', '') - if n_cyc: - cfg['index']['nterm_peptide_mods_spec'] = _apply_mod(spec, mod=NCYC_MOD) - else: - cfg['index']['nterm_peptide_mods_spec'] = _apply_mod(spec, mod='', exclusive_pattern=r'^\d*Q\-17\.027') - logMsg.debug(f'{'--n-cyc' if n_cyc else '--no-n-cyc'} applied: nterm_peptide_mods_spec updated to {cfg['index']['nterm_peptide_mods_spec']}') - if n_ace is not None: - spec = cfg['index'].get('nterm_protein_mods_spec', '') - if n_ace: - cfg['index']['nterm_protein_mods_spec'] = _apply_mod(spec, mod=NACE_MOD) - else: - cfg['index']['nterm_protein_mods_spec'] = _apply_mod(spec, mod='', exclusive_pattern=r'^\d*X\+42\.011') - logMsg.debug(f'{'--n-ace' if n_ace else '--no-n-ace'} applied: nterm_protein_mods_spec updated to {cfg['index']['nterm_protein_mods_spec']}') - if low_res is not None: - if low_res: - cfg['search']['mz_bin_width'] = MZ_BIN_WIDTH_LOW_RES - cfg['search']['score_function'] = SCORE_FUNC_LOW_RES - logMsg.debug(f'--low-res applied: mz_bin_width: {cfg["search"]["mz_bin_width"]}, score_function: {cfg["search"]["score_function"]}') - else: - cfg['search']['mz_bin_width'] = MZ_BIN_WIDTH_HIGH_RES - cfg['search']['score_function'] = SCORE_FUNC_HIGH_RES - logMsg.debug(f'--high-res applied: mz_bin_width: {cfg["search"]["mz_bin_width"]}, score_function: {cfg["search"]["score_function"]}') - logMsg.debug(f'{'--low-res' if low_res else '--high-res'} applied: mz_bin_width: {cfg["search"]["mz_bin_width"]}, score_function: {cfg["search"]["score_function"]}') - if clip_met is not None: - cfg.setdefault('index', {}) - cfg['index']['clip_n_met'] = 'true' if clip_met else 'false' - logMsg.debug(f'{'--clip-met' if clip_met else '--no-clip-met'} applied: {cfg['index']['clip_n_met']}') - return cfg - -# -- _apply_mod: returns mod_spec string -def _apply_mod(mods_spec: str, mod: str, exclusive_pattern: str | None = None) -> str: - ''' - Add or remove a mod entry in a Tide mods_spec string. - ''' - # Split on commas, discard empty strings from a blank mods_spec - entries = [e.strip() for e in mods_spec.split(',') if e.strip()] - if exclusive_pattern: - pattern = re.compile(exclusive_pattern, re.IGNORECASE) - entries = [e for e in entries if not pattern.match(e)] - elif mod == '': - pass - else: - entries = [e for e in entries if e != mod] - if mod: - entries = [mod] + entries - return ','.join(entries) - -# -- _apply-iodo: returns fixed_mods string -def _apply_iodo(fixed_mods: str, iodo: bool) -> str: - ''' - Add or remove the carbamidomethylation Cys mod in a Tide fixed_mods string - ''' - # Split on commas, discard empty strings from a blank mods_spec - entries = [e.strip() for e in fixed_mods.split(',') if e.strip()] - entries = [e for e in entries if e != CARBAMIDOMETHYL_MOD and e != 'C+0'] - if iodo: - entries = [CARBAMIDOMETHYL_MOD] + entries - else: - entries = ['C+0'] + entries # Crux automatically adds cysteine carbamidomethylation unless this string present - result = ','.join(entries) - return result - -def _apply_custom_mod(custom_mods: str, new_entry: str) -> str: - ''' - Add a custom mod entry to the custom_mods string, or clear all custom mods if new_entry is an empty string - ''' - if new_entry == '': - return '' - # Check against managed mod patterns - for pattern, flag_name in MANAGED_MOD_PATTERNS.items(): - if re.match(pattern, new_entry, re.IGNORECASE): - logMsg.warn(f'{new_entry} is managed by the {flag_name} flag, ignoring') - return custom_mods - # Split on commas, discard empty strings from a blank mods_spec - entries = [e.strip() for e in custom_mods.split(',') if e.strip()] - if new_entry not in entries: - entries.append(new_entry) - out_str = ','.join(entries) - logMsg.debug(f'custom_mods updated: {out_str}') - return out_str - -# -- _apply_organism: returns config dict with organism section replaced -def _apply_organism(cfg: dict, organism: dict[str, str]) -> dict: - ''' - Replace the [organism] section of the user config with the supplied dictionary. - ''' - cfg['organism'] = organism - logMsg.debug(f'organism section set to {organism}') - return cfg - -# -- _parse_organism_arg: returns dict parsed from list of 'Key=Pattern' strings -def _parse_organism_arg(pairs: list[str]) -> dict[str, str]: - ''' - Parse a list of 'Label=Pattern' strings into a dict. - ''' - result = {} - for item in pairs: - if '=' not in item: - logMsg.error(f'Invalid organism argument {item} (expected format: Organism=Pattern)') - raise SystemExit(1) - key, _, pattern = item.partition('=') - key = ''.join(key.split()) - pattern = ''.join(pattern.split()) - if not key: - logMsg.error(f'Empty label in organism argument: {item}') - raise SystemExit(1) - if not pattern: - logMsg.error(f'Empty pattern in organism argument: {item}') - raise SystemExit(1) - result[key] = pattern - return result - -# _mod_summary_line: prints a s -def _mod_summary_line(flag: bool | None, mod: str, key: str): - ''' - Print a single ✓ line for a boolean mod flag, or nothing if flag is None - ''' - if flag is None: - return - print(f'[bold green]✓[/bold green] [dim]{key}[/dim] → [cyan]{mod}[/cyan]') - -# _print_set_summary: prints a summary of changes made -def _printSetSummary( - *, - iodo: bool | None, - ox: bool | None, - phos: bool | None, - n_cyc: bool | None, - n_ace: bool | None, - low_res: bool | None, - organism: list[str] | None, - custom: str | None, - clip_met: bool | None, -) -> None: - ''' - Print a summary of what config_set changed - ''' - print() - _mod_summary_line(iodo, CARBAMIDOMETHYL_MOD, f'index.{'fixed_mods'}') - _mod_summary_line(ox, MET_OX_MOD, 'index.mods_spec') - _mod_summary_line(phos, PHOSPHO_MOD, 'index.mods_spec') - if custom is not None: - if custom == '': - print(f'[bold green]✓[/bold green] Custom mods cleared: [dim]index.custom_mods[/dim] → [cyan](empty)[/cyan]') - else: - print(f'[bold green]✓[/bold green] Custom mod added: [dim]index.custom_mods[/dim] → [cyan]{custom}[/cyan]') - _mod_summary_line(n_cyc, NCYC_MOD, f'index.{'nterm_peptide_mods_spec'}') - _mod_summary_line(n_ace, NACE_MOD, f'index.{'nterm_protein_mods_spec'}') - if clip_met is not None: - value = 'true' if clip_met else 'false' - print(f'[bold green]✓[/bold green] Clipped N-terminal methionine set: [dim]index.clip_n_met[/dim] → to [cyan]{value}[/cyan]') - if low_res is not None: - if low_res: - print(f'[bold green]✓[/bold green] Low-resolution mode set: [dim]search.mz_bin_width[/dim] → [cyan]{MZ_BIN_WIDTH_LOW_RES}[/cyan], [dim]search.score_function[/dim] → [cyan]{SCORE_FUNC_LOW_RES}[/cyan]') - else: - print(f'[bold green]✓[/bold green] High-resolution mode set: [dim]search.mz_bin_width[/dim] → [cyan]{MZ_BIN_WIDTH_HIGH_RES}[/cyan], [dim]search.score_function[/dim] → [cyan]{SCORE_FUNC_HIGH_RES}[/cyan]') - if organism is not None: - for item in organism: - key, _, pattern = item.partition('=') - key = ''.join(key.split()) - pattern = ''.join(pattern.split()) - print(f'[bold green]✓[/bold green] Organism pattern set: [dim]organism[/dim] → [cyan]{key}[/cyan]: [cyan]{pattern}[/cyan]') - print() \ No newline at end of file + _print_diff_summary(before, _flatten(cfg)) + return True \ No newline at end of file diff --git a/src/comms/commands/convert.py b/src/comms/commands/convert.py index fd32adf..51dae04 100644 --- a/src/comms/commands/convert.py +++ b/src/comms/commands/convert.py @@ -9,16 +9,22 @@ # -- Import internal functions from comms.utils.log import configureFileLogging, logMsg from comms.utils.context import ExperimentContext, resolve_data_files +from comms.utils.settings import resolve_config_value, _writeConfigTo from comms.utils.validate import validate from comms.utils import trfp as trfputil from comms.utils import paths as pathutil # -- run_convert: converts all .RAW files in input_dir to indexed mzML and writes them to output -def run_convert(data_files, ctx: ExperimentContext, gzip: bool | None = None, in_pipeline: bool = False): - if not in_pipeline: - logMsg('convert') +def run_convert( + data_files, + ctx: ExperimentContext, + gzip: bool | None = None, + format: int | None = None, + metadata: int | None = None, + in_pipeline: bool = False, +): + logMsg('convert') logMsg.debug('Started command: convert') - gzip = ctx.config['convert']['gzip'] if gzip is None else gzip _, trfp_path = validate(check_trfp=True, bin_dir=ctx.bin_dir) data_files = resolve_data_files(ctx, data_files) raw_files = [f for f in data_files if f.suffix.lower() == '.raw'] @@ -38,6 +44,19 @@ def run_convert(data_files, ctx: ExperimentContext, gzip: bool | None = None, in log_path = out_dir / 'convert.log' configureFileLogging(log_path) logMsg.debug(f'Output log file: {log_path}') + # Build config for this run only if any override was given + overrides_given = any(v is not None for v in (gzip, format, metadata)) + if overrides_given: + logMsg.debug('Using run-specific configuration parameters') + run_config = {**ctx.config, 'convert': dict(ctx.config.get('convert', {}))} + run_config['convert']['gzip'] = resolve_config_value(ctx.config, 'convert', 'gzip', gzip) + run_config['convert']['format'] = resolve_config_value(ctx.config, 'convert', 'format', format) + run_config['convert']['metadata'] = resolve_config_value(ctx.config, 'convert', 'metadata', metadata) + logMsg.info('Command-line overrides detected - run configuration file will be saved to output folder as "convert.config.toml"') + _writeConfigTo(run_config, path=Path(out_dir, 'convert.config.toml')) + else: + logMsg.debug('Using contextual configuration parameters') + run_config = ctx.config n_ok, n_fail = 0, 0 for raw_file in raw_files: logMsg.progress(f'Converting {raw_file.name}') @@ -45,8 +64,8 @@ def run_convert(data_files, ctx: ExperimentContext, gzip: bool | None = None, in trfp_path=trfp_path, raw_file=raw_file, out_dir=out_dir, - output_format=ctx.config['convert']['format'], - metadata=ctx.config['convert']['metadata'], + output_format=run_config['convert']['format'], + metadata=run_config['convert']['metadata'], ) if ok: n_ok += 1 @@ -71,7 +90,7 @@ def run_convert(data_files, ctx: ExperimentContext, gzip: bool | None = None, in except: continue # If --gzip was provided, gzip TRFP output - if gzip: + if run_config['convert']['gzip']: import gzip, os, shutil # Loop through each mzML file in directory for file in out_dir.glob('[!.]*.mzML'): diff --git a/src/comms/commands/experiment.py b/src/comms/commands/experiment.py index dbf1ff4..47d4c08 100644 --- a/src/comms/commands/experiment.py +++ b/src/comms/commands/experiment.py @@ -3,20 +3,40 @@ ''' # -- Import external dependencies -import tomli_w, typer +import re, tomli_w, tomllib from datetime import datetime, timezone from pathlib import Path from rich import print -from typing import Literal # -- Import internal functions +from comms.commands.config import _writeConfigTo +from comms.utils.context import _normalise_dirs +from comms.utils.installrdeps import check_r_dependencies, install_r_dependencies from comms.utils.log import logMsg -from comms.utils.sheet import SampleRow, render_sample_sheet -from comms.commands.config import _apply_protocol_flags, _apply_organism, _writeConfigTo +from comms.utils.modspec import apply_protocol_flags, apply_organism from comms.utils.settings import loadDefaultConfig +from comms.utils.sheet import SampleRow, render_sample_sheet, parse_sample_sheet + +# -- _existing_experiment: returns (root, comms_dir, metadata, config, sample_rows) if experiment_dir already holds a saved experiment else None +def _existing_experiment(experiment_dir: Path): + root, comms_dir = _normalise_dirs(experiment_dir) + root = Path(root).resolve() + meta_path = comms_dir / 'experiment.toml' + config_path = comms_dir / 'config.toml' + if not (meta_path.exists() and config_path.exists()): + return None + with meta_path.open('rb') as f: + metadata = tomllib.load(f) + with config_path.open('rb') as f: + config = tomllib.load(f) + sheet_path = comms_dir / 'sample_sheet.tsv' + rows: list[SampleRow] = [] + if sheet_path.exists(): + rows = parse_sample_sheet(sheet_path.read_text(encoding='utf-8')) + return root, comms_dir, metadata, config, rows # -- launch_experiment_gui: opens the PySide6 experiment setup window -def launch_experiment_gui() -> None: +def launch_experiment_gui(experiment_dir: Path | None = None) -> None: logMsg('experiment') try: from comms.gui.app import run_app @@ -24,13 +44,15 @@ def launch_experiment_gui() -> None: logMsg.error(f'Could not import GUI components: {e}') raise SystemExit(1) logMsg.info(f'Launching experiment setup GUI') - raise SystemExit(run_app()) + raise SystemExit(run_app(experiment_dir)) -# -- _prompt_list: return a list of strings by repeated prompting -def _prompt_list(label: str) -> list[str]: - items: list[str] = [] +# -- _prompt_list: return a list of strings by repeated prompting, prepopulated with any existing items +def _prompt_list(label: str, existing: list[str] | None = None) -> list[str]: + items: list[str] = list(existing or []) + if items: + print(f'Current {label}(s): {", ".join(items)}') while True: - value = typer.prompt(f'Add a {label} (blank to finish)', default='', show_default=False) + value = logMsg.input(f'Add a {label} (blank to finish)', default='', show_default=False) value = value.strip() if not value: break @@ -38,32 +60,49 @@ def _prompt_list(label: str) -> list[str]: items.append(value) return items -# -- _choose: prompt until the user picks one of the allowed options -def _choose(label: str, options: list[str]) -> str: - while True: - choice = typer.prompt(f'{label} {options}') - if choice in options: - return choice +# -- _choose: prompt until the user picks one of the allowed options, prepopulated with any existing values +def _choose(label: str, options: list[str], default: str | None = None) -> str: + formatted_options = [f'{opt} (default)' if opt==default else opt for opt in options] + label = f'{label} [dim]\\[{", ".join(formatted_options)}][/dim]' + return logMsg.input(label, choices=options, default=default, show_default=False, show_choices=False) -# -- run_experiment_headless: build a sample sheet, config and metadata via prompts -def run_experiment_headless() -> None: +# -- _confirm: yes/no prompt via logMsg.input, returned as a bool +def _confirm(msg: str, default: bool) -> bool: + msg = f'{msg} [dim]({"Y/n" if default else "y/N"})[/dim]' + answer = logMsg.input(msg, choices=['y', 'n'], default='y' if default else 'n', case_sensitive=False, show_choices=False, show_default=False) + return str(answer).strip().lower() == 'y' + +# -- run_experiment_headless: build a sample sheet, config and metadata via prompts, or edit an existing experiment +def run_experiment_headless(experiment_dir: Path | None = None) -> None: logMsg('experiment') logMsg.debug('Starting command: experiment') - logMsg.info('Starting headless experiment setup') - - name = typer.prompt('Experiment name') - base_dir = Path(typer.prompt('Save experiment to (directory)')).expanduser() - bin_dir = typer.prompt('Bin directory (blank to auto-resolve)', default='', show_default=False).strip() - database = typer.prompt('Combined database FASTA').strip() + existing = _existing_experiment(experiment_dir) if experiment_dir else None + edit_mode = existing is not None + if edit_mode: + root, comms_dir, metadata, config, existing_rows = existing + logMsg.info(f'Existing experiment found at {comms_dir}, editing in place') + else: + metadata, config, existing_rows = {}, {}, [] + logMsg.info('Starting headless experiment edit' if edit_mode else 'Starting headless experiment setup') - treatments = _prompt_list('treatment') - fractions = _prompt_list('fraction') + name = logMsg.input('Experiment name', default=metadata.get('experiment', {}).get('name', ''), show_default=edit_mode) + if edit_mode: + base_dir = root + else: + base_dir = Path(logMsg.input('Save experiment to directory (default: ".")', default=str(experiment_dir) if experiment_dir else '.', show_default=experiment_dir is not None)).expanduser() + bin_dir = logMsg.input('Bin directory (blank to auto-resolve)', default=metadata.get('experiment', {}).get('bin_dir', ''), show_default=edit_mode).strip() + database = logMsg.input('Combined database FASTA', default=metadata.get('files', {}).get('database', ''), show_default=edit_mode).strip() + existing_treatments = sorted({r.treatment for r in existing_rows if r.treatment}) + existing_fractions = sorted({r.fraction for r in existing_rows if r.fraction}) + treatments = _prompt_list('treatment', existing=existing_treatments) + fractions = _prompt_list('fraction', existing=existing_fractions) if not treatments or not fractions: logMsg.error('At least one treatment and one fraction are required') raise SystemExit(1) - input_dir = Path(typer.prompt('Directory of .RAW / .mzML files')).expanduser() - input_files = _prompt_list('data file') + input_dir = Path(logMsg.input('Directory of .RAW / .mzML files')).expanduser() + existing_data_files = metadata.get('files', {}).get('data', []) if edit_mode else None + input_files = _prompt_list('data file', existing=existing_data_files) files = [] for f in input_files: f = Path(Path(f).expanduser()) @@ -74,54 +113,67 @@ def run_experiment_headless() -> None: logMsg.error(f'No .RAW or .mzML files found in {input_dir}') raise SystemExit(1) + existing_by_raw = {r.raw_file: r for r in existing_rows} rows: list[SampleRow] = [] counters: dict[tuple[str, str], int] = {} for f in files: print(f'\n[bold]{f.name}[/bold]') - treatment = _choose('Treatment', treatments) - fraction = _choose('Fraction', fractions) + prior = existing_by_raw.get(f.name) + treatment = _choose('Treatment', treatments, default=prior.treatment if prior else None) + fraction = _choose('Fraction', fractions, default=prior.fraction if prior else None) key = (treatment, fraction) counters[key] = counters.get(key, 0) + 1 rows.append(SampleRow( - sample_id=f.stem, raw_file=f.name, - treatment=treatment, fraction=fraction, replicate=counters[key], + sample_id=prior.sample_id if prior else f.stem, + raw_file=f.name, + treatment=treatment, + fraction=fraction, + replicate=counters[key], )) # Config: reuse the same helpers as the GUI's ConfigPanel + index_cfg = config.get('index', {}) + search_cfg = config.get('search', {}) cfg = loadDefaultConfig() - cfg = _apply_protocol_flags( + cfg = apply_protocol_flags( cfg, - iodo=typer.confirm('Cysteine carbamidomethylation (static)?', default=False), - ox=typer.confirm('Methionine oxidation (variable)?', default=True), - phos=typer.confirm('STY phosphorylation (variable)?', default=False), - n_cyc=typer.confirm('N-terminal Gln cyclisation?', default=True), - n_ace=typer.confirm('Protein N-terminal acetylation?', default=True), - clip_met=typer.confirm('Clip N-terminal methionine?', default=True), - low_res=typer.confirm('Low-resolution instrument (ion trap)?', default=False), + iodo=_confirm('Cysteine carbamidomethylation (static)?', default='C+0' not in index_cfg.get('fixed_mods', '')), + ox=_confirm('Methionine oxidation (variable)?', default=bool(re.search(r'M\+15\.9949', index_cfg.get('mods_spec', ''))) if edit_mode else True), + phos=_confirm('STY phosphorylation (variable)?', default=bool(re.search(r'STY\+79\.966331', index_cfg.get('mods_spec', ''))) if edit_mode else False), + n_cyc=_confirm('N-terminal Gln cyclisation?', default=bool(index_cfg.get('nterm_peptide_mods_spec', '')) if edit_mode else True), + n_ace=_confirm('Protein N-terminal acetylation?', default=bool(index_cfg.get('nterm_protein_mods_spec', '')) if edit_mode else True), + clip_met=_confirm('Clip N-terminal methionine?', default=(str(index_cfg.get('clip_n_met', True)).strip().lower() == 'true') if edit_mode else True), + low_res=_confirm('Low-resolution instrument (ion trap)?', default=(search_cfg.get('score_function') == 'combined-p-value') if edit_mode else False), ) - cfg.setdefault('index', {})['custom_mods'] = '' - # Analysis mode: single- or multi-species - organisms: dict[str, str] = {} - multispecies = typer.confirm('Multispecies analysis (per-organism FDR)?', default=False) + cfg.setdefault('index', {})['custom_mods'] = index_cfg.get('custom_mods', '') + organisms: dict[str, str] = dict(config.get('organism', {})) if edit_mode else {} + multispecies = _confirm('Multispecies analysis (per-organism FDR)?', default=bool(organisms)) if multispecies: while True: - label = typer.prompt('Organism label (blank to finish)', default='', show_default=False).strip() + label = logMsg.input('Organism label (blank to finish)', default='', show_default=False).strip() if not label: break - pattern = typer.prompt(f'Header pattern for {label}').strip() + pattern = logMsg.input(f'Header pattern for {label}').strip() if pattern: organisms[label] = pattern - cfg = _apply_organism(cfg, organisms) + else: + organisms = {} + cfg = apply_organism(cfg, organisms) if multispecies: - cfg['percolator']['shared_psm'] = typer.prompt(f'Shared PSM handling policy', default='drop', type=Literal['drop','include'], show_choices=True, show_default=True).strip() - # Report settings + cfg['rescore']['shared_psm'] = logMsg.input('Shared PSM handling policy', choices=['drop', 'include'], default=config.get('rescore', {}).get('shared_psm', 'drop'), show_choices=True, show_default=True).strip() + report_meta = metadata.get('report', {}) organism_prefix = '' - include_report = typer.confirm('Create report?', default=True) + include_report = _confirm('Create report?', default=report_meta.get('enabled', True)) if include_report: - reference = typer.prompt('Reference protein annotation file (blank to skip)', default='', show_default=False).strip() - contaminants = typer.prompt('Contaminants list CSV path (blank to skip)', default='', show_default=False).strip() + reference = logMsg.input('Reference protein annotation file (blank to skip)', default=report_meta.get('ref_info', ''), show_default=edit_mode).strip() + contaminants = logMsg.input('Contaminants list CSV path (blank to skip)', default=report_meta.get('cont_csv', ''), show_default=edit_mode).strip() if multispecies: - organism_prefix = typer.prompt('Primary organism ID prefix', default='', show_default=False).strip() + organism_prefix = logMsg.input('Primary organism ID prefix', default=report_meta.get('organism_prefix', ''), show_default=edit_mode).strip() + status = check_r_dependencies() + if status is not None and status['missing']: + if _confirm(f"Install missing R report dependencies now? ({', '.join(status['missing'])})", default=True): + install_r_dependencies() + # Write all three files out_dir = base_dir / 'comms' out_dir.mkdir(parents=True, exist_ok=True) @@ -153,9 +205,8 @@ def run_experiment_headless() -> None: meta['report']['organism_prefix'] = organism_prefix else: meta.setdefault('report', {})['enabled'] = False - # Save metadata file with (out_dir / 'experiment.toml').open('wb') as f: tomli_w.dump(meta, f) - logMsg.info(f'Experiment written to {out_dir}') + logMsg.info(f'Experiment {"updated" if edit_mode else "written"} at {out_dir}') print(f'\nRun the pipeline with:\n' f'\t[bold]comms pipeline {sheet_path} --database --input {input_dir} --experiment-dir {base_dir}[/bold]\n') \ No newline at end of file diff --git a/src/comms/commands/index.py b/src/comms/commands/index.py index 6ce5446..685df72 100644 --- a/src/comms/commands/index.py +++ b/src/comms/commands/index.py @@ -9,14 +9,26 @@ # -- Import internal functions from comms.utils.log import configureFileLogging, logMsg from comms.utils.context import ExperimentContext, resolve_database +from comms.utils.modspec import apply_protocol_flags, apply_custom_mod +from comms.utils.settings import _writeConfigTo from comms.utils.validate import validate from comms.utils import crux as cruxutil from comms.utils import paths as pathutil -# -- run_index: builds a Tide peptide index from database and writes it to output -def run_index(database, ctx: ExperimentContext, in_pipeline: bool = False): - if not in_pipeline: - logMsg('index') +def run_index( + database, + ctx: ExperimentContext, + in_pipeline: bool = False, + iodo=None, + ox=None, + phos=None, + n_cyc=None, + n_ace=None, + custom=None, + clip_met=None, + missed_cleavages=None, +): + logMsg('index') logMsg.debug('Started command: index') crux_bin, _ = validate(check_crux=True, bin_dir=ctx.bin_dir) database = resolve_database(ctx, database) @@ -26,12 +38,35 @@ def run_index(database, ctx: ExperimentContext, in_pipeline: bool = False): log_path = out_dir / 'index.log' configureFileLogging(log_path) logMsg.debug(f'Output log file: {log_path}') + # Build a config for this run only if any override was given + overrides_given = any(v is not None for v in (iodo, ox, phos, n_cyc, n_ace, custom, clip_met, missed_cleavages)) + if overrides_given: + logMsg.debug('Using run-specific configuration parameters') + run_config = {**ctx.config, 'index': dict(ctx.config.get('index', {}))} + run_config = apply_protocol_flags( + run_config, + iodo=iodo, + ox=ox, + phos=phos, + n_cyc=n_cyc, + n_ace=n_ace, + clip_met=clip_met, + missed_cleavages=missed_cleavages, + ) + if custom is not None: + current = run_config['index'].get('custom_mods', '') + run_config['index']['custom_mods'] = apply_custom_mod(current, custom) + logMsg.info('Command-line overrides detected - run configuration file will be saved to output folder as "index.config.toml"') + _writeConfigTo(run_config, path=Path(out_dir, 'index.config.toml')) + else: + logMsg.debug('Using contextual configuration parameters') + run_config = ctx.config logMsg.progress(f'Building Tide peptide index') ok = cruxutil.tideIndex( crux_bin=crux_bin, database=database, index_dir=out_dir, - config=ctx.config, + config=run_config, ) if not ok: logMsg.error(f'tide-index failed, see {log_path}') diff --git a/src/comms/commands/lfq.py b/src/comms/commands/lfq.py index 31e4013..c41fb66 100644 --- a/src/comms/commands/lfq.py +++ b/src/comms/commands/lfq.py @@ -23,8 +23,7 @@ def run_lfq( ctx: ExperimentContext, in_pipeline: bool = False ): - if not in_pipeline: - logMsg('lfq') + logMsg('lfq') logMsg.debug('Started command: lfq') crux_bin, _ = validate(check_crux=True, allow_lfq=True, bin_dir=ctx.bin_dir) rescore_dir = resolve_results_input(ctx, 'rescore', rescore_dir) @@ -79,4 +78,4 @@ def _groupPsmsByFraction(psm_files: list[Path], samples: pd.DataFrame) -> dict[s # -- _get_stem: return str corresponding to stem of file from raw_file in sample sheet def _get_stem(row): - return str(row['sample_id']) \ No newline at end of file + return Path(str(row['raw_file'])).stem \ No newline at end of file diff --git a/src/comms/commands/pipeline.py b/src/comms/commands/pipeline.py index 965ad06..4c1491a 100644 --- a/src/comms/commands/pipeline.py +++ b/src/comms/commands/pipeline.py @@ -7,7 +7,7 @@ from rich import print # -- Import internal functions -from comms.utils.log import logMsg +from comms.utils.log import logMsg, concatenatePipelineLog, startPipelineLogging from comms.utils.samples import loadSampleSheet from comms.utils.context import ExperimentContext, resolve_database, resolve_data_files, resolve_report, resolve_sample_sheet from comms.commands import convert, index, search, rescore, lfq, quantify, report @@ -28,6 +28,7 @@ def run_pipeline( ): logMsg('pipeline') logMsg.debug(f'Started command: pipeline') + startPipelineLogging() # Resolve external inputs once data_files = resolve_data_files(ctx, data) @@ -47,19 +48,27 @@ def run_pipeline( raise SystemExit(1) logMsg.debug(f"Sample sheet loaded: {len(samples)} sample(s); {samples['treatment'].nunique()} treatment(s)") logMsg.info(f"Running comMS pipeline: {len(samples)} sample(s), {samples['treatment'].nunique()} treatment(s)") + # -- Step 1: Convert (optional) if not skip_convert: current_step += 1 logMsg.progress(f'Step {current_step}/{num_steps}: converting .RAW files') - convert.run_convert(data_files, ctx=ctx, gzip=True, in_pipeline=True) + convert.run_convert( + data_files, + ctx=ctx, + gzip=None, + in_pipeline=True + ) mzml_override = None # search/lfq glob the convert results else: logMsg.progress(f'Skipped .RAW -> .mzML conversion') mzml_override = [f for f in data_files if f.suffix.lower() == '.mzml' or f.name.endswith('.mzML.gz')] + # -- Step 2: Build index current_step += 1 logMsg.progress(f'Step {current_step}/{num_steps}: building peptide index') index.run_index(database=database, ctx=ctx, in_pipeline=True) + # -- Step 3: Search current_step += 1 logMsg.progress(f'Step {current_step}/{num_steps}: searching spectra') @@ -71,6 +80,7 @@ def run_pipeline( threads=threads, in_pipeline=True ) + # -- Step 4: Rescore current_step += 1 logMsg.progress(f'Step {current_step}/{num_steps}: rescoring PSMs') @@ -80,7 +90,8 @@ def run_pipeline( ctx=ctx, organism_tags=org_tags, in_pipeline=True -) + ) + # -- Steps 5 & 6: Quantify if skip_lfq and skip_quantify: logMsg.progress(f'Skipped LFQ and dNSAF quantification') @@ -112,17 +123,21 @@ def run_pipeline( ref_info=None, cont_csv=None, organism_prefix=None, - # ! TODO: make below configurable via CLI or config? - min_reps=3, - fdr_threshold=0.05, - lfc_threshold=1.0, + min_reps=None, + fdr_threshold=None, + lfc_threshold=None, + top_n=None, sections=VALID_SECTIONS, overwrite=False, rscript='Rscript', ) END = datetime.datetime.now() + logMsg('pipeline') # logger needs to be re-tagged logMsg.info(f'Pipeline complete, runtime {END - START}, results written to {ctx.root}') + pipeline_log_path = concatenatePipelineLog(ctx.root / 'comms/results/pipeline.log') + if pipeline_log_path is not None: + logMsg.debug(f'Aggregated pipeline log written to: {pipeline_log_path}') logMsg.debug(f'Finished command: pipeline') # -- _calculate_n_steps: returns int corresponding to number of steps in pipeline diff --git a/src/comms/commands/quantify.py b/src/comms/commands/quantify.py index 75cb4c3..43dfa10 100644 --- a/src/comms/commands/quantify.py +++ b/src/comms/commands/quantify.py @@ -11,14 +11,22 @@ # -- Import internal functions from comms.utils.log import configureFileLogging, logMsg from comms.utils.context import ExperimentContext, resolve_database, resolve_results_input +from comms.utils.settings import resolve_config_value, _writeConfigTo from comms.utils.validate import validate from comms.utils import crux as cruxutil from comms.utils import paths as pathutil -# -- run_quantify: runs dNSAF spectral counting on rescored PSM files (discovering single/multi-species results) and writes results to output -def run_quantify(input_dir, database, ctx: ExperimentContext, in_pipeline: bool = False): - if not in_pipeline: - logMsg('quantify') +# -- run_quantify: run spectral counting on rescored PSM files (discovering single/multi-species results) and writes results to output +def run_quantify( + input_dir, + database, + ctx: ExperimentContext, + measure=None, + qvalue_threshold=None, + unique_mapping=None, + in_pipeline: bool = False, +): + logMsg('quantify') logMsg.debug('Started command: quantify') crux_bin, _ = validate(check_crux=True, bin_dir=ctx.bin_dir) input_dir = resolve_results_input(ctx, 'rescore', input_dir) @@ -40,6 +48,21 @@ def run_quantify(input_dir, database, ctx: ExperimentContext, in_pipeline: bool log_path = out_dir / 'quantify.log' configureFileLogging(log_path) logMsg.debug(f'Output log file: {log_path}') + + # Build config for this run only if any override was given + overrides_given = any(v is not None for v in (measure, qvalue_threshold, unique_mapping)) + if overrides_given: + logMsg.debug('Using run-specific configuration parameters') + run_config = {**ctx.config, 'quantify': dict(ctx.config.get('quantify', {}))} + run_config['quantify']['measure'] = resolve_config_value(ctx.config, 'quantify', 'measure', measure) + run_config['quantify']['qvalue_threshold'] = resolve_config_value(ctx.config, 'quantify', 'qvalue_threshold', qvalue_threshold) + run_config['quantify']['unique_mapping'] = resolve_config_value(ctx.config, 'quantify', 'unique_mapping', unique_mapping) + logMsg.info('Command-line overrides detected - run configuration file will be saved to output folder as "quantify.config.toml"') + _writeConfigTo(run_config, path=Path(out_dir, 'quantify.config.toml')) + else: + logMsg.debug('Using contextual configuration parameters') + run_config = ctx.config + n_ok, n_fail = 0, 0 with logging_redirect_tqdm(): for psm_file in tqdm(psm_files, desc='Files quantified'): @@ -51,7 +74,7 @@ def run_quantify(input_dir, database, ctx: ExperimentContext, in_pipeline: bool database=database, out_dir=out_dir, fileroot=fileroot, - config=ctx.config, + config=run_config, ) if ok: n_ok += 1 diff --git a/src/comms/commands/report.py b/src/comms/commands/report.py index 05b2952..b75eb29 100644 --- a/src/comms/commands/report.py +++ b/src/comms/commands/report.py @@ -3,7 +3,7 @@ ''' # -- Import external dependencies -import shutil, subprocess, sys +import json, shutil, subprocess, sys from datetime import datetime from importlib.resources import files as pkg_files from pathlib import Path @@ -11,31 +11,72 @@ from rich.console import Console # -- Import internal functions -from comms.utils.log import logMsg +from comms.utils.log import configureFileLogging, logMsg from comms.utils.samples import loadSampleSheet +from comms.utils.settings import resolve_config_value, _writeConfigTo from comms.utils.context import ExperimentContext, resolve_organism_prefix, resolve_sample_sheet, resolve_results_input, results_dir # -- Initialise Rich console console = Console() -# -- Define helper dictionary matching sections to R scripts and whether they require LFQ data -_SECTIONS: dict[str, tuple[str, bool]] = { - # Core sections - 'qc': ('qc.R', False), - 'pca': ('pca.R', False), - 'da': ('da.R', False), - 'secondary-species': ('secondary-species.R', False), - 'concordance': ('concordance.R', True), +# -- Define helper dictionary matching sections to R scripts, LFQ requirement, and whether the section reports per-organism status +_SECTIONS: dict[str, tuple[str, bool, bool]] = { + # Core sections: (script, needs_lfq, per_organism) + 'qc': ('qc.R', False, True), + 'pca': ('pca.R', False, True), + 'da': ('da.R', False, True), + 'secondary-species': ('secondary-species.R', False, False), + 'concordance': ('concordance.R', True, True), # Auxiliary sections - 'ev-markers': ('aux/ev-markers.R', False), + 'ev-markers': ('aux/ev-markers.R', False, True), } +# -- Define helper dictionary for potential per-organism section statuses +_ORGANISM_STATUSES = {'ok', 'skipped', 'failed'} + # -- _resolve_r_script: returns Path to R script def _resolve_r_script(script_name: str) -> Path: '''Locate R script based on provided script name''' return pkg_files('comms').joinpath(f'r/sections/{script_name}') -# -- _run_r_script: returns boolean indicating if command was run successfully +# -- _read_status: returns tuple of dicts (containing organisms, reasons) parsed from a section's _status.json +def _read_status(output_subdir: Path) -> tuple[dict[str, str], dict[str, str]]: + status_path = output_subdir / '_status.json' + if not status_path.exists(): + return {}, {} + try: + payload = json.loads(status_path.read_text()) + except Exception as e: + logMsg.debug(f'Could not parse {status_path}: {e}') + return {}, {} + organisms = {k: v for k, v in payload.get('organisms', {}).items() if v in _ORGANISM_STATUSES} + reasons = dict(payload.get('reasons', {})) + return organisms, reasons + +# -- _section_status: returns string (either 'succeeded', 'partial', 'failed' or 'skipped') +def _section_status(proc_ok: bool, organisms: dict[str, str]) -> str: + if not organisms: + # No structured status available (legacy/non-organism script, or crash before writing status) + return 'failed' if not proc_ok else 'skipped' + ok = sum(1 for s in organisms.values() if s == 'ok') + failed = sum(1 for s in organisms.values() if s == 'failed') + if failed == 0: + return 'succeeded' if ok > 0 else 'skipped' + return 'partial' if ok > 0 else 'failed' + +# -- _log_organism_outcomes: returns None but outputs logging messages +def _log_organism_outcomes(section: str, organisms: dict[str, str], reasons: dict[str, str]) -> None: + for org, status in organisms.items(): + reason = reasons.get(org) + suffix = f' ({reason})' if reason else '' + if status == 'ok': + logMsg.info(f'{section} — {org}: succeeded') + elif status == 'skipped': + logMsg.info(f'{section} — {org}: skipped{suffix}') + else: + logMsg.warn(f'{section} — {org}: failed{suffix}') + +# -- _run_r_section: returns boolean indicating if the R process itself exited cleanly def _run_r_section( section: str, script_name: str, @@ -56,7 +97,12 @@ def _run_r_section( return True # -- _write_index: return none, but write index -def _write_index(output_dir: Path, params: dict, results: dict[str, bool]) -> None: +def _write_index( + output_dir: Path, + params: dict, + section_status: dict[str, str], + organism_results: dict[str, dict[str, str]], +) -> None: lines = [ '# comms report', f'\nGenerated: {datetime.now().strftime("%Y-%m-%d %H:%M:%S")}', @@ -65,9 +111,16 @@ def _write_index(output_dir: Path, params: dict, results: dict[str, bool]) -> No for k,v in params.items(): lines.append(f'- **{k}**: `{v}`') lines.append(f'\n## Sections\n') - for sec, ok in results.items(): - status = '✓ SUCCEEDED' if ok else '✗ FAILED' - lines.append(f'- {sec}: {status}') + status_glyphs = { + 'succeeded': '✓ SUCCEEDED', + 'partial': '◐ PARTIAL', + 'failed': '✗ FAILED', + 'skipped': '- SKIPPED', + } + for sec, status in section_status.items(): + lines.append(f'- {sec}: {status_glyphs[status]}') + for org, org_status in organism_results.get(sec, {}).items(): + lines.append(f' - {org}: {org_status}') (output_dir / 'index.md').write_text('\n'.join(lines)) # -- run_report: return None, but run report section R scripts and output script @@ -79,19 +132,22 @@ def run_report( ref_info: Path | None, cont_csv: Path | None, organism_prefix: str | None, - min_reps: int, - lfc_threshold: float, - fdr_threshold: float, + min_reps: int | None, + lfc_threshold: float | None, + fdr_threshold: float | None, + top_n: int | None, sections: list, overwrite: bool, rscript: str, in_pipeline: bool, ) -> None: - if not in_pipeline: - logMsg('report') + logMsg('report') logMsg.debug('Started command: report') # Create path to output_dir output_dir = ctx.root / 'comms/results/report' + log_path = output_dir / 'report.log' + configureFileLogging(log_path) + logMsg.debug(f'Output log file: {log_path}') # Required inputs resolved from the experiment context quantify_dir = resolve_results_input(ctx, 'quantify', quantify_dir) @@ -136,6 +192,22 @@ def run_report( if shutil.which(rscript) is None: logMsg.error(f'Rscript not callable: {rscript}') raise SystemExit(1) + + # Build config for this run only if any override was given + overrides_given = any(v is not None for v in (min_reps, lfc_threshold, fdr_threshold, top_n)) + if overrides_given: + logMsg.debug('Using run-specific configuration parameters') + run_config = {**ctx.config, 'report': dict(ctx.config.get('report', {}))} + run_config['report']['min_reps'] = resolve_config_value(ctx.config, 'report', 'min_reps', min_reps) + run_config['report']['lfc_threshold'] = resolve_config_value(ctx.config, 'report', 'lfc_threshold', lfc_threshold) + run_config['report']['fdr_threshold'] = resolve_config_value(ctx.config, 'report', 'fdr_threshold', fdr_threshold) + run_config['report']['top_n_proteins'] = resolve_config_value(ctx.config, 'report', 'top_n_proteins', top_n) + logMsg.info('Command-line overrides detected - run configuration file will be saved to output folder as "report.config.toml"') + _writeConfigTo(run_config, path=Path(output_dir, 'report.config.toml')) + else: + logMsg.debug('Using contextual configuration parameters') + run_config = ctx.config + # Run command logMsg.info(f'Generating report: {len(sections)} section(s)') # Define arguments passed to every R script @@ -145,24 +217,32 @@ def run_report( str(ref_info) if ref_info else '', str(cont_csv) if cont_csv else '', organism_prefix, - str(min_reps), + str(run_config['report']['min_reps']), ] - results: dict[str, bool] = {} + section_status: dict[str, str] = {} + organism_results: dict[str, dict[str, str]] = {} for sec in sections: logMsg.progress(f'Running section: {sec}') - script, needs_lfq = _SECTIONS[sec] + script, needs_lfq, per_organism = _SECTIONS[sec] extra: list[str] = [] if sec == 'da': - extra = [str(lfc_threshold), str(fdr_threshold)] + extra = [str(run_config['report']['lfc_threshold']), str(run_config['report']['fdr_threshold']), str(run_config['report']['top_n_proteins'])] elif sec == 'concordance': - extra = [str(lfq_dir), str(lfc_threshold), str(fdr_threshold)] - results[sec] = _run_r_section( + extra = [str(lfq_dir), str(run_config['report']['lfc_threshold']), str(run_config['report']['fdr_threshold'])] + output_subdir = output_dir / sec.replace('-', '_') + proc_ok = _run_r_section( section = sec, script_name=script, - output_subdir=output_dir / sec.replace('-', '_'), + output_subdir=output_subdir, positional_args=common_args+extra, rscript=rscript, ) + organisms, reasons = _read_status(output_subdir) if per_organism else ({}, {}) + section_status[sec] = _section_status(proc_ok, organisms) + organism_results[sec] = organisms + if organisms: + _log_organism_outcomes(sec, organisms, reasons) + _write_index( output_dir, { @@ -170,12 +250,25 @@ def run_report( 'sample_sheet': sample_sheet, 'lqf_dir': lfq_dir or 'not provided', 'organism_prefix': organism_prefix, - 'min_reps': min_reps, - 'lfc_threshold': lfc_threshold, - 'fdr_threshold': fdr_threshold, + 'min_reps': run_config['report']['min_reps'], + 'lfc_threshold': run_config['report']['lfc_threshold'], + 'fdr_threshold': run_config['report']['fdr_threshold'], + 'top_n_proteins': run_config['report']['top_n_proteins'], }, - results) - n_ok = sum(results.values()) - n_fail = len(results) - n_ok - logMsg.info(f'Report complete: {n_ok} succeeded, {n_fail} failed') + section_status, + organism_results, + ) + + n_succeeded = sum(1 for s in section_status.values() if s == 'succeeded') + n_partial = sum(1 for s in section_status.values() if s == 'partial') + n_failed = sum(1 for s in section_status.values() if s == 'failed') + n_skipped = sum(1 for s in section_status.values() if s == 'skipped') + parts = [f'{n_succeeded} succeeded'] + if n_partial: + parts.append(f'{n_partial} partial (at least one organism failed)') + if n_failed: + parts.append(f'{n_failed} failed') + if n_skipped: + parts.append(f'{n_skipped} skipped (no organism had sufficient data)') + logMsg.info(f'Report complete: {", ".join(parts)}') logMsg.debug(f'Finished command: report') \ No newline at end of file diff --git a/src/comms/commands/rescore.py b/src/comms/commands/rescore.py index aa7fd54..3c7f8c5 100644 --- a/src/comms/commands/rescore.py +++ b/src/comms/commands/rescore.py @@ -13,6 +13,7 @@ from comms.utils.fasta import splitFastaByOrganism from comms.utils.log import configureFileLogging, logMsg from comms.utils.context import ExperimentContext, resolve_database, resolve_results_input +from comms.utils.settings import resolve_config_value, _writeConfigTo from comms.utils.validate import validate from comms.utils import crux as cruxutil from comms.utils import paths as pathutil @@ -25,10 +26,12 @@ def run_rescore( database, ctx: ExperimentContext, organism_tags: Optional[str] = None, + protein_enzyme: Optional[str] = None, + picked_protein: Optional[bool] = None, + shared_psm: Optional[str] = None, in_pipeline: bool = False, ): - if not in_pipeline: - logMsg('rescore') + logMsg('rescore') logMsg.debug('Started command: rescore') crux_bin, _ = validate(check_crux=True, bin_dir=ctx.bin_dir) input_dir = resolve_results_input(ctx, 'search', input_dir) @@ -51,20 +54,36 @@ def run_rescore( log_path = out_dir / 'rescore.log' configureFileLogging(log_path) logMsg.debug(f'Output log file: {log_path}') + + # Build config for this run only if any override was given + overrides_given = any(v is not None for v in (protein_enzyme, picked_protein, shared_psm)) + if overrides_given: + logMsg.debug('Using run-specific configuration parameters') + run_config = {**ctx.config, 'rescore': dict(ctx.config.get('rescore', {}))} + run_config['rescore']['protein_enzyme'] = resolve_config_value(ctx.config, 'rescore', 'protein_enzyme', protein_enzyme) + run_config['rescore']['picked_protein'] = resolve_config_value(ctx.config, 'rescore', 'picked_protein', picked_protein) + run_config['rescore']['shared_psm'] = resolve_config_value(ctx.config, 'rescore', 'shared_psm', shared_psm) + logMsg.info('Command-line overrides detected - run configuration file will be saved to output folder as "rescore.config.toml"') + _writeConfigTo(run_config, path=Path(out_dir, 'rescore.config.toml')) + else: + logMsg.debug('Using contextual configuration parameters') + run_config = ctx.config + # Round 1: run Percolator on the full combined database, with one call per sample file combined_target_files = _run_combined_percolator_round( - crux_bin, target_files, database, out_dir, ctx, + crux_bin, + target_files, + database, + out_dir, + run_config, ) if not combined_target_files: logMsg.error('No combined Percolator output found, cannot continue') raise SystemExit(1) + # Round 2: run Percolator on each organism sub-FASTA if multispecies analysis if multispecies: - organism_tags = ( - _parseOrganismTags(organism_tags) - if organism_tags - else ctx.config.get('organism') - ) + organism_tags = (_parseOrganismTags(organism_tags) if organism_tags else ctx.config.get('organism')) if not organism_tags: logMsg.error('No organism tags supplied or configured for multi-species analysis') raise SystemExit(1) @@ -72,13 +91,18 @@ def run_rescore( sub_fastas = splitFastaByOrganism(database, out_dir, organism_tags) logMsg.debug(f'Built {len(sub_fastas)} per-organism sub-FASTA(s)') _run_per_organism_percolator_round( - crux_bin, target_files, sub_fastas, organism_tags, out_dir, ctx + crux_bin, + target_files, + sub_fastas, + organism_tags, + out_dir, + run_config, ) # Log command as complete logMsg.debug('Finished command: rescore') # -- _run_percolator_round: returns list of PSM files after runn Percolator (via Crux) on database, with one call per sample file -def _run_combined_percolator_round(crux_bin, target_files, database, out_dir, ctx) -> list: +def _run_combined_percolator_round(crux_bin, target_files, database, out_dir, run_config) -> list: logMsg.progress(f'Rescoring {len(target_files)} file(s) using combined database') n_ok, n_fail = 0, 0 with logging_redirect_tqdm(): @@ -91,7 +115,7 @@ def _run_combined_percolator_round(crux_bin, target_files, database, out_dir, ct database=database, out_dir=out_dir, fileroot=fileroot, - config=ctx.config, + config=run_config, ) if ok: n_ok += 1 @@ -103,10 +127,10 @@ def _run_combined_percolator_round(crux_bin, target_files, database, out_dir, ct return sorted(out_dir.glob('[!.]*.percolator.target.psms.txt')) # -- _run_per_organism_percolator_round: returns None but splits combined Tide search outputs by organism and runs Percolator (via Crux) -def _run_per_organism_percolator_round(crux_bin, combined_target_files, sub_fastas, organism_tags, out_dir, ctx): +def _run_per_organism_percolator_round(crux_bin, combined_target_files, sub_fastas, organism_tags, out_dir, run_config): logMsg.progress(f'Rescoring {len(combined_target_files)} file(s) using per-organism sub-FASTAs') n_ok, n_fail = 0, 0 - shared_policy = ctx.config['percolator']['shared_psm'] + shared_policy = run_config['rescore']['shared_psm'] with logging_redirect_tqdm(): for combined_file in tqdm(combined_target_files, desc='Files rescored'): logMsg.progress(f'Rescoring {combined_file.name}') @@ -140,7 +164,7 @@ def _run_per_organism_percolator_round(crux_bin, combined_target_files, sub_fast database=sub_fastas[label], out_dir=org_out_dir, fileroot=org_fileroot, - config=ctx.config, + config=run_config, ) if ok: n_ok += 1 diff --git a/src/comms/commands/search.py b/src/comms/commands/search.py index 6cb634e..e65060a 100644 --- a/src/comms/commands/search.py +++ b/src/comms/commands/search.py @@ -11,34 +11,61 @@ # -- Import internal functions from comms.utils.log import configureFileLogging, logMsg from comms.utils.context import ExperimentContext, resolve_results_input, resolve_mzml_files +from comms.utils.settings import resolve_config_value, _writeConfigTo from comms.utils.validate import validate from comms.utils import crux as cruxutil from comms.utils import paths as pathutil # -- run_search: runs tide-search on all mzML files in input_dir and writes results to output -def run_search(data_files, index_dir, ctx: ExperimentContext, param_medic: bool, threads: int, in_pipeline: bool = False): - if not in_pipeline: - logMsg('search') +def run_search( + data_files, + index_dir, + ctx: ExperimentContext, + param_medic: bool, + threads: int | None, + score_function: str | None = None, + min_peaks: int | None = None, + precursor_tolerance_ppm: float | None = None, + mz_bin_width: float | None = None, + in_pipeline: bool = False, +): + logMsg('search') logMsg.debug('Started command: search') crux_bin, _ = validate(check_crux=True, bin_dir=ctx.bin_dir) index_dir = resolve_results_input(ctx, 'index', index_dir) mzml_files = resolve_mzml_files(ctx, data_files) - threads = threads or ctx.config['search']['threads'] logMsg.info(f'Searching {len(mzml_files)} mzML file(s)') out_dir = pathutil.generateOutputFileStructure(ctx.root, 'search') logMsg.debug(f'Output directory: {out_dir}') log_path = out_dir / 'search.log' configureFileLogging(log_path) logMsg.debug(f'Output log file: {log_path}') + # -- Optional: param-medic tolerance estimation - precursor_tol = None - mz_bin_width = None + pm_precursor, pm_bin_width = None, None if param_medic: logMsg.progress(f'Estimating tolerances with param-medic') - precursor_tol, mz_bin_width = _runParamMedic(crux_bin=crux_bin, mzml_files=mzml_files, out_dir=out_dir) - prec_display = precursor_tol or ctx.config['search']['precursor_tolerance_ppm'] - bin_width_display = mz_bin_width or ctx.config['search']['mz_bin_width'] - logMsg.debug(f'Precursor tolerance {prec_display} ppm, m/z bin width {bin_width_display} Da') + pm_precursor, pm_bin_width = _runParamMedic(crux_bin=crux_bin, mzml_files=mzml_files, out_dir=out_dir) + + # Build config for this run only if any override was given + overrides_given = any(v is not None for v in (threads, score_function, min_peaks, precursor_tolerance_ppm, mz_bin_width)) + if overrides_given: + logMsg.debug('Using run-specific configuration parameters') + run_config = {**ctx.config, 'search': dict(ctx.config.get('search', {}))} + run_config['search']['threads'] = resolve_config_value(ctx.config, 'search', 'threads', threads) + run_config['search']['score_function'] = resolve_config_value(ctx.config, 'search', 'score_function', score_function) + run_config['search']['min_peaks'] = resolve_config_value(ctx.config, 'search', 'min_peaks', min_peaks) + run_config['search']['precursor_tolerance_ppm'] = resolve_config_value(ctx.config, 'search', 'precursor_tolerance_ppm', precursor_tolerance_ppm if precursor_tolerance_ppm is not None else pm_precursor) + run_config['search']['mz_bin_width'] = resolve_config_value(ctx.config, 'search', 'mz_bin_width', mz_bin_width if mz_bin_width is not None else pm_bin_width) + logMsg.info('Command-line overrides detected - run configuration file will be saved to output folder as "search.config.toml"') + _writeConfigTo(run_config, path=Path(out_dir, 'search.config.toml')) + else: + logMsg.debug('Using contextual configuration parameters') + run_config = ctx.config + + logMsg.debug(f'Precursor tolerance {run_config["search"]["precursor_tolerance_ppm"]} ppm, m/z bin width {run_config["search"]["mz_bin_width"]} Da') + + # -- PSM search n_ok, n_fail = 0, 0 with logging_redirect_tqdm(): for mzml_file in tqdm(mzml_files, desc='Files searched'): @@ -50,10 +77,7 @@ def run_search(data_files, index_dir, ctx: ExperimentContext, param_medic: bool, index_dir=index_dir, out_dir=out_dir, fileroot=fileroot, - config=ctx.config, - threads=threads, - precursor_tol=prec_display, - mz_bin_width=bin_width_display, + config=run_config, ) if ok: n_ok += 1 diff --git a/src/comms/config.toml b/src/comms/config.toml index bff566f..c9b0ba3 100644 --- a/src/comms/config.toml +++ b/src/comms/config.toml @@ -1,9 +1,5 @@ # comMS default configuration -[global] -verbose = false -debug = false - [organism] @@ -22,30 +18,24 @@ nterm_protein_mods_spec = "1X+42.011" custom_mods = "" [search] -score_function = "xcorr" # xcorr for high-res MS (default); combined-p-value for low-res (run: comms config set --low-res) -mz_bin_width = 0.02 # 0.02 for high-res MS (default); 1.0005079 for low-res (run: comms config set --low-res) +score_function = "xcorr" # xcorr for high-res MS (default); combined-p-value for low-res (run: comms config --low-res) +mz_bin_width = 0.02 # 0.02 for high-res MS (default); 1.0005079 for low-res (run: comms config --low-res) min_peaks = 10 precursor_tolerance_ppm = 10.0 threads = 2 -[percolator] +[rescore] protein_enzyme = "trypsin" -picked_protein = true # use picked-protein FDR (c.f. Savitski et al. 2015) +picked_protein = true # use picked-protein FDR (c.f. Savitski et al. 2015, doi:10.1074/mcp.M114.046995) shared_psm = "drop" -[lfq] -match_between_runs = true - [quantify] measure = "dNSAF" qvalue_threshold = 0.01 unique_mapping = true [report] -top_n_proteins = 50 -fdr_threshold = 0.01 -colour_palette = "Set2" -include_organism_panel = true -include_pca = true -include_volcano = true -include_heatmap = true \ No newline at end of file +min_reps = 3 +lfc_threshold = 1.0 +fdr_threshold = 0.05 +top_n_proteins = 20 \ No newline at end of file diff --git a/src/comms/gui/app.py b/src/comms/gui/app.py index d4bf967..8bda882 100644 --- a/src/comms/gui/app.py +++ b/src/comms/gui/app.py @@ -4,14 +4,15 @@ # -- Import external dependencies import sys +from pathlib import Path from PySide6.QtWidgets import QApplication # -- Import internal functions from comms.gui.main_window import MainWindow # -- run_app: create the QApplication, show the main window, and run the event loop -def run_app() -> int: +def run_app(experiment_dir: Path | None = None) -> int: app = QApplication.instance() or QApplication(sys.argv) - window = MainWindow() + window = MainWindow(experiment_dir=experiment_dir) window.show() return app.exec() \ No newline at end of file diff --git a/src/comms/gui/main_window.py b/src/comms/gui/main_window.py index 3d3d97d..d35bcf6 100644 --- a/src/comms/gui/main_window.py +++ b/src/comms/gui/main_window.py @@ -3,6 +3,7 @@ ''' # -- Import external dependencies +from pathlib import Path from PySide6.QtCore import QSize from PySide6.QtWidgets import QMainWindow, QTabWidget, QVBoxLayout, QWidget @@ -18,7 +19,7 @@ # -- MainWindow: four numbered tabs with per-tab status icons class MainWindow(QMainWindow): - def __init__(self, parent=None): + def __init__(self, experiment_dir: Path | None = None, parent=None): super().__init__(parent) self._log = logMsg('experiment') self.setWindowTitle('comms experiment setup') @@ -62,17 +63,39 @@ def __init__(self, parent=None): self.sample.contentChanged.connect(self.readiness.refresh) self.config.changed.connect(self.readiness.refresh) self.experiment.changed.connect(self.readiness.refresh) + self.experiment.binDirChanged.connect(self.readiness.refresh_dependencies) # paint the initial (unedited) icons self.tabs.setTabIcon(self._sample_index, status_icon(self.sample.tracker.status)) self.tabs.setTabIcon(self._config_index, status_icon(self.config.tracker.status)) self.tabs.setTabIcon(self._review_index, status_icon(self.experiment.tracker.status)) self.save.refresh() + + # if a directory provided, call _load_existing + if experiment_dir is not None: + self._load_existing(experiment_dir) def _on_tab_changed(self, index: int) -> None: if index == self._review_index: self.save.refresh() - + + def _load_existing(self, experiment_dir: Path) -> None: + from comms.commands.experiment import _existing_experiment + existing = _existing_experiment(experiment_dir) + if existing is None: + # If no experiment actually saved to directory, only pre-fill save to field + self.experiment._dir.setText(str(experiment_dir)) + return + root, comms_dir, metadata, config, rows = existing + self.experiment.load_from_metadata(root, metadata) + report_meta = metadata.get('report', {}) + self.config.load_from_config(config, report_meta) + treatments = sorted({r.treatment for r in rows if r.treatment}) + fractions = sorted({r.fraction for r in rows if r.fraction}) + data_files = metadata.get('files', {}).get('data', []) + self.sample.load(rows, treatments, fractions, data_files=data_files) + self._log.info(f'Loaded existing experiment from {comms_dir}') + def closeEvent(self, event) -> None: self._log.info('Closed experiment setup GUI') super().closeEvent(event) \ No newline at end of file diff --git a/src/comms/gui/panels/config_panel.py b/src/comms/gui/panels/config_panel.py index 07f6327..78428ea 100644 --- a/src/comms/gui/panels/config_panel.py +++ b/src/comms/gui/panels/config_panel.py @@ -3,18 +3,20 @@ ''' # -- Import external dependencies +import re from pathlib import Path from PySide6.QtCore import Qt, Signal from PySide6.QtWidgets import ( QWidget, QVBoxLayout, QHBoxLayout, QFileDialog, QFormLayout, QGroupBox, QCheckBox, - QComboBox, QLineEdit, QPushButton, QTableWidget, + QComboBox, QLineEdit, QPushButton, QTableWidget, QTableWidgetItem ) # -- Import internal functions +from comms.commands.config import _writeConfigTo from comms.gui.status import PanelStateTracker from comms.utils.settings import loadDefaultConfig -from comms.commands.config import ( - _apply_protocol_flags, _apply_organism, _apply_custom_mod, _writeConfigTo, +from comms.utils.modspec import ( + apply_protocol_flags, apply_organism, apply_custom_mod, MZ_BIN_WIDTH_LOW_RES, ) # -- Define class ConfigPanel to define a structured form mirroring `comms config set` with an additional analysis type activating the organism table @@ -287,7 +289,7 @@ def summary(self) -> str: # -- build config file and save -- def _build_config(self) -> dict: cfg = loadDefaultConfig() - cfg = _apply_protocol_flags( + cfg = apply_protocol_flags( cfg, iodo=self._iodo.isChecked(), ox=self._ox.isChecked(), @@ -302,14 +304,14 @@ def _build_config(self) -> dict: organisms = { label: pattern for label, pattern in self._organism_rows() if label and pattern } - cfg['percolator']['shared_psm'] = self.shared_policy() - cfg = _apply_organism(cfg, organisms) + cfg['rescore']['shared_psm'] = self.shared_policy() + cfg = apply_organism(cfg, organisms) cfg.setdefault('index', {}) cfg['index']['custom_mods'] = '' custom = self._custom.text().strip() if custom: for entry in [e.strip() for e in custom.split(',') if e.strip()]: - cfg['index']['custom_mods'] = _apply_custom_mod( + cfg['index']['custom_mods'] = apply_custom_mod( cfg['index']['custom_mods'], entry) return cfg @@ -320,4 +322,38 @@ def sync_tracker(self) -> None: def write(self, out_dir: Path) -> Path: path = out_dir / 'config.toml' _writeConfigTo(self._build_config(), path) - return path \ No newline at end of file + return path + + # -- load_from_config: populate fields from a loaded config.toml + experiment.toml [report] section + def load_from_config(self, cfg: dict, report_meta: dict) -> None: + index_cfg = cfg.get('index', {}) + search_cfg = cfg.get('search', {}) + self._iodo.setChecked('C+0' not in index_cfg.get('fixed_mods', '')) + self._ox.setChecked(bool(re.search(r'M\+15\.9949', index_cfg.get('mods_spec', '')))) + self._phos.setChecked(bool(re.search(r'STY\+79\.966331', index_cfg.get('mods_spec', '')))) + self._n_cyc.setChecked(bool(index_cfg.get('nterm_peptide_mods_spec', ''))) + self._n_ace.setChecked(bool(index_cfg.get('nterm_protein_mods_spec', ''))) + clip_met_value = index_cfg.get('clip_n_met', True) + if isinstance(clip_met_value, str): + clip_met_value = clip_met_value.strip().lower() == 'true' + self._clip_met.setChecked(bool(clip_met_value)) + self._custom.setText(index_cfg.get('custom_mods', '')) + self._res.setCurrentIndex(1 if search_cfg.get('mz_bin_width') == MZ_BIN_WIDTH_LOW_RES else 0) + + organisms = cfg.get('organism', {}) + self._analysis.setCurrentIndex(1 if organisms else 0) + self._sharedpsm.setCurrentIndex(1 if cfg.get('rescore', {}).get('shared_psm') == 'include' else 0) + self._org_table.setRowCount(0) + for label, pattern in organisms.items(): + row = self._org_table.rowCount() + self._org_table.insertRow(row) + self._org_table.setItem(row, 0, QTableWidgetItem(label)) + self._org_table.setItem(row, 1, QTableWidgetItem(pattern)) + + self._report_enabled.setChecked(bool(report_meta.get('enabled', True))) + self._reference.setText(str(report_meta.get('ref_info', ''))) + self._contaminant.setText(str(report_meta.get('cont_csv', ''))) + self._organism_prefix.setText(str(report_meta.get('organism_prefix', ''))) + self._update_organism_enabled() + self._update_report_fields_enabled() + self._on_changed() \ No newline at end of file diff --git a/src/comms/gui/panels/experiment_panel.py b/src/comms/gui/panels/experiment_panel.py index 937cd76..b01e054 100644 --- a/src/comms/gui/panels/experiment_panel.py +++ b/src/comms/gui/panels/experiment_panel.py @@ -17,6 +17,7 @@ # -- Define class ExperimentPanel to collect experiment name and base output directory class ExperimentPanel(QWidget): changed = Signal() + binDirChanged = Signal() def __init__(self, parent=None): super().__init__(parent) @@ -59,6 +60,7 @@ def __init__(self, parent=None): self._bin.setMinimumWidth(360) self._bin.setPlaceholderText('optional: directory containing Crux / ThermoRawFileParser') self._bin.textChanged.connect(self.changed) + self._bin.editingFinished.connect(self.binDirChanged) bin_browse = QPushButton('Select directory') bin_browse.clicked.connect(self._browse_bin) bin_row = QWidget() @@ -117,6 +119,12 @@ def bin_dir(self) -> Path | None: text = self._bin.text().strip() return Path(text) if text else None + # -- set_bin_dir: write a chosen bin directory into the field and mark the panel changed + def set_bin_dir(self, path: Path) -> None: + self._bin.setText(str(path)) + self.changed.emit() + self.binDirChanged.emit() + def database_path(self) -> Path | None: text = self._database.text().strip() return Path(text) if text else None @@ -148,5 +156,19 @@ def write_metadata(self, out_dir: Path, files: dict | None = None, analysis=None tomli_w.dump(meta, f) return path + # -- load_from_metadata: populate fields from a loaded experiment.toml + resolved base_dir + def load_from_metadata(self, base_dir: Path, metadata: dict) -> None: + self._name.setText(metadata.get('experiment', {}).get('name', '')) + self._dir.setText(str(base_dir)) + database = metadata.get('files', {}).get('database', '') + if database: + self._database.setText(str(database)) + bin_dir = metadata.get('experiment', {}).get('bin_dir', '') + if bin_dir: + self._bin.setText(str(bin_dir)) + self.changed.emit() + if bin_dir: + self.binDirChanged.emit() + def is_valid(self) -> bool: return bool(self.experiment_name()) and self.base_dir() is not None and self.database_path() is not None \ No newline at end of file diff --git a/src/comms/gui/panels/readiness_panel.py b/src/comms/gui/panels/readiness_panel.py index 8c33c2c..90be63e 100644 --- a/src/comms/gui/panels/readiness_panel.py +++ b/src/comms/gui/panels/readiness_panel.py @@ -3,13 +3,17 @@ ''' # -- Import external dependencies +from pathlib import Path from PySide6.QtCore import Qt from PySide6.QtWidgets import ( - QWidget, QVBoxLayout, QGridLayout, QGroupBox, QLabel, + QApplication, QFileDialog, QGridLayout, QGroupBox, QHBoxLayout, QLabel, QMessageBox, QPushButton, QVBoxLayout, QWidget, ) # -- Import internal functions +from comms.utils.installrdeps import check_r_dependencies, install_r_dependencies +from comms.utils.paths import repoBinDir from comms.utils.readiness import COMMANDS, missing_requirements +from comms.utils.validate import probe_crux, probe_trfp from comms.gui.widgets.status_indicator import StatusIndicator from comms.gui.status import PanelStatus @@ -24,7 +28,10 @@ def __init__(self, experiment, sample, config, parent=None): self._experiment = experiment self._sample = sample self._config = config - + # -- cached dependency state, refreshed by _refresh_dependencies() rather than on every refresh() + self._r_deps_status: dict | None = None + self._crux_found = False + self._trfp_found = False layout = QVBoxLayout(self) box = QGroupBox('Command readiness') grid = QGridLayout(box) @@ -43,11 +50,34 @@ def __init__(self, experiment, sample, config, parent=None): self._indicators[command] = indicator self._details[command] = detail + deps_box = QGroupBox('Dependencies') + deps_layout = QVBoxLayout(deps_box) + summary_row = QHBoxLayout() + self._deps_indicator = StatusIndicator() + self._deps_label = QLabel() + summary_row.addWidget(self._deps_indicator, 0, Qt.AlignmentFlag.AlignVCenter) + summary_row.addWidget(self._deps_label, 1) + deps_layout.addLayout(summary_row) + + buttons_row = QHBoxLayout() + self._install_deps_btn = QPushButton('Install R dependencies') + self._install_deps_btn.clicked.connect(self._install_deps) + self._locate_trfp_btn = QPushButton('Locate TRFP…') + self._locate_trfp_btn.clicked.connect(self._set_bin_dir) + self._locate_crux_btn = QPushButton('Locate Crux…') + self._locate_crux_btn.clicked.connect(self._set_bin_dir) + buttons_row.addWidget(self._install_deps_btn) + buttons_row.addWidget(self._locate_trfp_btn) + buttons_row.addWidget(self._locate_crux_btn) + deps_layout.addLayout(buttons_row) + + layout.addWidget(deps_box) layout.addWidget(box) layout.addStretch(1) + self._refresh_dependencies() self.refresh() - # -- _state: gather the booleans the readiness model needs from the source panels + # -- _state: gather the booleans the readiness model needs from the source panels and cached dependency state def _state(self) -> dict: return dict( has_data=len(self._sample.data_files()) > 0, @@ -56,13 +86,16 @@ def _state(self) -> dict: has_organism_prefix=bool(self._config.organism_prefix()), multispecies=self._config.analysis_mode() == 'multi', has_organism_tags=self._config.has_organism_patterns(), + has_trfp=self._trfp_found, + has_crux=self._crux_found, + has_r_deps=self._r_deps_status is not None and not self._r_deps_status['missing'], ) - # -- refresh: recompute readiness and repaint each row + # -- refresh: recompute command-row readiness and repaint each row (reads cached dependency state; call refresh_dependencies() to re-probe) def refresh(self) -> None: missing = missing_requirements(**self._state()) for command in COMMANDS: - gaps = list(dict.fromkeys(missing[command])) # dedupe, keep order (pipeline repeats) + gaps = list(dict.fromkeys(missing[command])) ready = not gaps self._indicators[command].setStatus( PanelStatus.COMPLETE if ready else PanelStatus.UNEDITED @@ -71,4 +104,67 @@ def refresh(self) -> None: self._details[command].setText(text) tooltip = 'Ready to run' if ready else 'Missing: ' + ', '.join(gaps) self._indicators[command].setToolTip(tooltip) - self._details[command].setToolTip(tooltip) \ No newline at end of file + self._details[command].setToolTip(tooltip) + + # -- refresh_dependencies: public entry point for main_window.py to call after the bin_dir field changes + def refresh_dependencies(self) -> None: + self._refresh_dependencies() + self.refresh() + + # -- _refresh_dependencies: re-probe R deps / Crux / TRFP, repaint the Dependencies box, and cache results for _state() + def _refresh_dependencies(self) -> None: + self._r_deps_status = check_r_dependencies() + bin_dir = repoBinDir(experiment_bin_dir=self._experiment.bin_dir()) + self._crux_found = probe_crux(bin_dir) is not None + self._trfp_found = probe_trfp(bin_dir) is not None + parts = [] + if self._r_deps_status is None: + parts.append('Rscript not found on PATH') + elif self._r_deps_status['missing']: + parts.append(f"{len(self._r_deps_status['missing'])} R package(s) missing") + else: + parts.append('R packages OK') + parts.append('Crux OK' if self._crux_found else 'Crux not found') + parts.append('TRFP OK' if self._trfp_found else 'TRFP not found') + self._deps_label.setText(' · '.join(parts)) + tooltip_lines = [f'Searched for Crux/TRFP in: {bin_dir}'] + if self._r_deps_status and self._r_deps_status['missing']: + tooltip_lines.append(f"Missing R packages: {', '.join(self._r_deps_status['missing'])}") + tooltip = '\n'.join(tooltip_lines) + self._deps_label.setToolTip(tooltip) + r_ok = self._r_deps_status is not None and not self._r_deps_status['missing'] + all_ok = r_ok and self._crux_found and self._trfp_found + none_ok = not r_ok and not self._crux_found and not self._trfp_found + if all_ok: + self._deps_indicator.setStatus(PanelStatus.COMPLETE) + elif none_ok: + self._deps_indicator.setStatus(PanelStatus.UNEDITED) + else: + self._deps_indicator.setStatus(PanelStatus.INCOMPLETE) + self._deps_indicator.setToolTip(tooltip) + + can_install = self._r_deps_status is not None and bool(self._r_deps_status['missing']) + self._install_deps_btn.setEnabled(can_install) + + # -- _set_bin_dir: open a directory picker and write the chosen path into the experiment panel's bin_dir field + def _set_bin_dir(self) -> None: + current = self._experiment.bin_dir() + start_dir = str(current) if current else '' + chosen = QFileDialog.getExistingDirectory( + self, 'Select bin directory (containing Crux and/or ThermoRawFileParser)', start_dir, + ) + if not chosen: + return + self._experiment.set_bin_dir(Path(chosen)) + + def _install_deps(self) -> None: + QApplication.setOverrideCursor(Qt.CursorShape.WaitCursor) + try: + ok = install_r_dependencies() + finally: + QApplication.restoreOverrideCursor() + self.refresh_dependencies() + if ok: + QMessageBox.information(self, 'R dependencies', 'All comms report R dependencies are installed.') + else: + QMessageBox.warning(self, 'R dependencies', 'Some R dependencies could not be installed. Check the terminal log for details.') \ No newline at end of file diff --git a/src/comms/gui/panels/sample_panel.py b/src/comms/gui/panels/sample_panel.py index 957442e..77b7573 100644 --- a/src/comms/gui/panels/sample_panel.py +++ b/src/comms/gui/panels/sample_panel.py @@ -58,4 +58,19 @@ def sync_tracker(self) -> None: def write(self, out_dir: Path) -> Path: path = out_dir / 'sample_sheet.tsv' path.write_text(self.sample_sheet_text(), encoding='utf-8') - return path \ No newline at end of file + return path + + # -- load: populate the sample table and treatment/fraction groups from a loaded sample sheet + def load(self, rows: list, treatments: list[str], fractions: list[str], data_files: list[str] | None = None) -> None: + for t in treatments: + self._state.add_treatment(t) + for f in fractions: + self._state.add_fraction(f) + if data_files: + by_name = {Path(p).name: str(p) for p in data_files} + for row in rows: + if not row.source_path: + match = by_name.get(row.raw_file) + if match: + row.source_path = match + self._state.sample_model.set_rows(rows) \ No newline at end of file diff --git a/src/comms/r/deps/check_deps.R b/src/comms/r/deps/check_deps.R new file mode 100644 index 0000000..7a77f70 --- /dev/null +++ b/src/comms/r/deps/check_deps.R @@ -0,0 +1,22 @@ +#!/bin/R +# check_deps.R: report which R dependencies are installed (as JSON to stdout) + +script_dir <- local({ + args <- commandArgs(trailingOnly = FALSE) + script <- grep("^--file=", args, value = TRUE) + dirname(normalizePath(sub("^--file=", "", script))) +}) +source(file.path(script_dir, "dependencies.R")) + +all_packages <- c(R_DEPENDENCIES$cran, R_DEPENDENCIES$bioc) +installed <- character(0) +missing <- character(0) +for (pkg in all_packages) { + ok <- tryCatch(requireNamespace(pkg, quietly = TRUE), error = function(e) FALSE) + if (ok) installed <- c(installed, pkg) else missing <- c(missing, pkg) +} +cat(sprintf( + '{"installed":[%s],"missing":[%s]}', + paste(sprintf('"%s"', installed), collapse = ","), + paste(sprintf('"%s"', missing), collapse = ",") +)) \ No newline at end of file diff --git a/src/comms/r/deps/dependencies.R b/src/comms/r/deps/dependencies.R new file mode 100644 index 0000000..a6f0cc3 --- /dev/null +++ b/src/comms/r/deps/dependencies.R @@ -0,0 +1,7 @@ +#!/bin/R +# dependencies.R: list of all R packages required by comms report command + +R_DEPENDENCIES <- list( + cran = c("tidyverse", "openxlsx2", "svglite", "ggrepel", "ggfortify", "cluster", "UpSetR", "pheatmap", "VennDiagram", "iq", "jsonlite"), + bioc = c("limma") +) \ No newline at end of file diff --git a/src/comms/r/deps/install_deps.R b/src/comms/r/deps/install_deps.R new file mode 100644 index 0000000..1deaed9 --- /dev/null +++ b/src/comms/r/deps/install_deps.R @@ -0,0 +1,37 @@ +#!/bin/R +# install_deps.R: install R dependencies for comms report command (only installing missing packages) + +script_dir <- local({ + args <- commandArgs(trailingOnly = FALSE) + script <- grep("^--file=", args, value = TRUE) + dirname(normalizePath(sub("^--file=", "", script))) +}) +source(file.path(script_dir, "dependencies.R")) + +is_missing <- function(pkg) !requireNamespace(pkg, quietly = TRUE) + +missing_cran <- Filter(is_missing, R_DEPENDENCIES$cran) +missing_bioc <- Filter(is_missing, R_DEPENDENCIES$bioc) + +if (length(missing_cran) == 0 && length(missing_bioc) == 0) { + message("All R report dependencies are already installed.") + quit(status = 0) +} + +if (length(missing_cran) > 0) { + message(sprintf("Installing %d CRAN package(s): %s", length(missing_cran), paste(missing_cran, collapse = ", "))) + install.packages(missing_cran, repos = "https://cloud.r-project.org") +} +if (length(missing_bioc) > 0) { + message(sprintf("Installing %d Bioconductor package(s): %s", length(missing_bioc), paste(missing_bioc, collapse = ", "))) + if (!requireNamespace("BiocManager", quietly = TRUE)) install.packages("BiocManager", repos = "https://cloud.r-project.org") + BiocManager::install(missing_bioc, ask = FALSE, update = FALSE) +} + +still_missing <- Filter(is_missing, c(R_DEPENDENCIES$cran, R_DEPENDENCIES$bioc)) +if (length(still_missing) > 0) { + message(sprintf("Still missing after install attempt: %s", paste(still_missing, collapse = ", "))) + quit(status = 1) +} + +message("All R report dependencies installed successfully.") \ No newline at end of file diff --git a/src/comms/r/install_deps.R b/src/comms/r/install_deps.R deleted file mode 100644 index 492ed15..0000000 --- a/src/comms/r/install_deps.R +++ /dev/null @@ -1,9 +0,0 @@ -#!/bin/R -# install_deps.R: install all required R dependencies for comms report command - -cran_packages <- c("tidyverse", "openxlsx2", "svglite", "ggrepel", "ggfortify", "cluster", "UpSetR", "pheatmap", "VennDiagram", "iq") -bioc_packages <- c("limma") - -install.packages(cran_packages, repos = "https://cloud.r-project.org") -if (!requireNamespace("BiocManager", quietly = TRUE)) install.packages("BiocManager", repos = "https://cloud.r-project.org") -BiocManager::install(bioc_packages, ask = FALSE, update = FALSE) \ No newline at end of file diff --git a/src/comms/r/sections/aux/ev-markers.R b/src/comms/r/sections/aux/ev-markers.R index b2d7c82..0afb7be 100644 --- a/src/comms/r/sections/aux/ev-markers.R +++ b/src/comms/r/sections/aux/ev-markers.R @@ -21,6 +21,7 @@ script_dir <- local({ # Import utility functions source(file.path(script_dir, "..", "..", "utils", "import.R")) source(file.path(script_dir, "..", "..", "utils", "normalise.R")) +source(file.path(script_dir, "..", "..", "utils", "state.R")) source(file.path(script_dir, "..", "..", "utils", "theme.R")) # Load libraries @@ -56,94 +57,107 @@ categorise_marker <- function(annotation) { organisms <- unique(sample_meta$organism) all_marker_tables <- list() +status <- new_status_tracker("ev-markers") for (org in organisms) { - org_meta <- filter(sample_meta, organism == org) - org_cols <- org_meta$dnsaf_col - # EV markers are only meaningful for the primary organism: skip any organism whose dNSAF columns contain no primary-organism proteins - primary_check <- results_wide %>% - filter(startsWith(proteinId, organism_prefix)) %>% - filter(rowSums(select(., all_of(org_cols)) > 0) > 0) - - if (nrow(primary_check) == 0) { - message(sprintf("EV markers: no primary organism proteins in %s columns — skipping", org)) - next - } - - fractions <- unique(org_meta$fraction) - - org_data <- results_wide %>% - filter(startsWith(proteinId, organism_prefix)) %>% - filter(rowSums(select(., all_of(org_cols)) > 0) > 0) - - # Per-fraction mean dNSAF - for (frac in fractions) { - cols <- intersect(filter(org_meta, fraction == frac)$dnsaf_col, colnames(org_data)) - org_data[[paste0("avg_", frac)]] <- - if (length(cols) > 0) rowMeans(select(org_data, all_of(cols)), na.rm=TRUE) else NA_real_ - } - - marker_table <- org_data %>% - rowwise() %>% - mutate(MISEVCategory=categorise_marker(proteinAnnotation)) %>% - ungroup() %>% - filter(!is.na(MISEVCategory)) %>% - mutate(MISEVCategory=factor(MISEVCategory, levels=MISEV_LEVELS)) %>% - arrange(MISEVCategory) - - if (nrow(marker_table) == 0) { - message(sprintf("EV markers %s: no marker proteins found — skipping", org)); next - } - - # Enrichment ratios - ev_frac <- fractions[str_detect(tolower(fractions), "ev")][1] - wcl_frac <- fractions[str_detect(tolower(fractions), "wcl")][1] - awf_frac <- fractions[str_detect(tolower(fractions), "awf|cr")][1] - if (!is.na(ev_frac) && !is.na(wcl_frac)) - marker_table <- mutate(marker_table, log2_EV_vs_WCL=log2((.data[[paste0("avg_", ev_frac)]] + 1e-10) / (.data[[paste0("avg_", wcl_frac)]] + 1e-10)) - if (!is.na(ev_frac) && !is.na(awf_frac)) - marker_table <- mutate(marker_table, log2_EV_vs_AWF=log2((.data[[paste0("avg_", ev_frac)]] + 1e-10) / (.data[[paste0("avg_", awf_frac)]] + 1e-10))) - - # Per-protein heatmap with category gaps - avg_cols <- intersect(paste0("avg_", fractions), colnames(marker_table)) - heatmap_mat <- marker_table %>% - select(proteinAnnotation, all_of(avg_cols)) %>% - column_to_rownames("proteinId") %>% - as.matrix() %>% - logdNSAF() - colnames(heatmap_mat) <- str_remove(colnames(heatmap_mat), "avg_") - ann_row <- data.frame(Category=as.character(marker_table$MISEVCategory), row.names=marker_table$proteinAnnotation) - gaps_row <- marker_table %>% - count(MISEVCategory) %>% - arrange(MISEVCategory) %>% - pull(n) %>% - cumsum() %>% - head(-1) - svglite(file.path(output_dir, sprintf("marker_heatmap_%s.svg", org)), width=12, height=max(6, nrow(heatmap_mat) * 0.35)) - pheatmap(heatmap_mat, annotation_row=ann_row, gaps_row=gaps_row, cluster_rows=FALSE, colour=colorRampPalette(c("#88CCEE", "white", "#CC6677"))(50), main=sprintf("MISEV2023 markers — log(dNSAF) — %s", org), fontsize_row=8, fontsize_col=10, border_colour=NA) - dev.off() - - # Aggregated category heatmap (3 × fraction×treatment) - agg_mat <- marker_table %>% - select(MISEVCategory, all_of(org_cols)) %>% - pivot_longer(-MISEVCategory, names_to="dnsaf_col", values_to="dNSAF") %>% - left_join(select(org_meta, dnsaf_col, fraction, treatment), by="dnsaf_col") %>% - mutate(log_dNSAF=logdNSAF(dNSAF)) %>% - group_by(MISEVCategory, fraction, treatment) %>% - summarise(mean_log_dNSAF=mean(log_dNSAF, na.rm=TRUE), .groups="drop") %>% - mutate(col_label=paste(fraction, treatment, sep="_")) %>% - select(MISEVCategory, col_label, mean_log_dNSAF) %>% - pivot_wider(names_from=col_label, values_from=mean_log_dNSAF) %>% - arrange(MISEVCategory) %>% - column_to_rownames("MISEVCategory") %>% - as.matrix() - - svglite(file.path(output_dir, sprintf("marker_category_heatmap_%s.svg", org)), width=8, height=4) - pheatmap(agg_mat, cluster_rows=FALSE, cluster_cols=FALSE, colour=colorRampPalette(c("#88CCEE", "white", "#CC6677"))(50), main=sprintf("Mean log(dNSAF) by MISEV category — %s", org), fontsize=10, border_colour=NA) - dev.off() - - all_marker_tables[[org]] <- marker_table + tryCatch({ + org_meta <- filter(sample_meta, organism == org) + org_cols <- org_meta$dnsaf_col + + primary_check <- results_wide %>% + filter(startsWith(proteinId, organism_prefix)) %>% + filter(rowSums(select(., all_of(org_cols)) > 0) > 0) + + if (nrow(primary_check) == 0) { + reason <- "no primary-organism proteins detected in this organism's samples" + message(sprintf("EV markers %s: %s — skipping", org, reason)) + status <<- record_skip(status, org, reason) + next + } + + fractions <- unique(org_meta$fraction) + org_data <- results_wide %>% + filter(startsWith(proteinId, organism_prefix)) %>% + filter(rowSums(select(., all_of(org_cols)) > 0) > 0) + + for (frac in fractions) { + cols <- intersect(filter(org_meta, fraction == frac)$dnsaf_col, colnames(org_data)) + org_data[[paste0("avg_", frac)]] <- + if (length(cols) > 0) rowMeans(select(org_data, all_of(cols)), na.rm = TRUE) else NA_real_ + } + + marker_table <- org_data %>% + rowwise() %>% + mutate(MISEVCategory=categorise_marker(proteinAnnotation)) %>% + ungroup() %>% + filter(!is.na(MISEVCategory)) %>% + mutate(MISEVCategory=factor(MISEVCategory, levels=MISEV_LEVELS)) %>% + arrange(MISEVCategory) + + if (nrow(marker_table) == 0) { + reason <- "no marker proteins found" + message(sprintf("EV markers %s: %s — skipping", org, reason)) + status <<- record_skip(status, org, reason) + next + } + + ev_frac <- fractions[str_detect(tolower(fractions), "ev")][1] + wcl_frac <- fractions[str_detect(tolower(fractions), "wcl")][1] + awf_frac <- fractions[str_detect(tolower(fractions), "awf|cr")][1] + if (!is.na(ev_frac) && !is.na(wcl_frac)) + marker_table <- mutate(marker_table, log2_EV_vs_WCL=log2((.data[[paste0("avg_", ev_frac)]] + 1e-10) / (.data[[paste0("avg_", wcl_frac)]] + 1e-10))) + if (!is.na(ev_frac) && !is.na(awf_frac)) + marker_table <- mutate(marker_table, log2_EV_vs_AWF=log2((.data[[paste0("avg_", ev_frac)]] + 1e-10) / (.data[[paste0("avg_", awf_frac)]] + 1e-10))) + + avg_cols <- intersect(paste0("avg_", fractions), colnames(marker_table)) + heatmap_mat <- marker_table %>% + select(proteinAnnotation, all_of(avg_cols)) %>% + column_to_rownames("proteinId") %>% + as.matrix() %>% + logdNSAF() + colnames(heatmap_mat) <- str_remove(colnames(heatmap_mat), "avg_") + + ann_row <- data.frame(Category=as.character(marker_table$MISEVCategory), row.names=marker_table$proteinAnnotation) + gaps_row <- marker_table %>% + count(MISEVCategory) %>% + arrange(MISEVCategory) %>% + pull(n) %>% + cumsum() %>% + head(-1) + + svglite(file.path(output_dir, sprintf("marker_heatmap_%s.svg", org)), width=12, height=max(6, nrow(heatmap_mat) * 0.35)) + pheatmap(heatmap_mat, annotation_row=ann_row, gaps_row=gaps_row, cluster_rows=FALSE, colour=colorRampPalette(c("#88CCEE", "white", "#CC6677"))(50), main=sprintf("MISEV2023 markers — log(dNSAF) — %s", org), fontsize_row=8, fontsize_col=10, border_colour=NA) + dev.off() + + agg_mat <- marker_table %>% + select(MISEVCategory, all_of(org_cols)) %>% + pivot_longer(-MISEVCategory, names_to="dnsaf_col", values_to="dNSAF") %>% + left_join(select(org_meta, dnsaf_col, fraction, treatment), by="dnsaf_col") %>% + mutate(log_dNSAF=logdNSAF(dNSAF)) %>% + group_by(MISEVCategory, fraction, treatment) %>% + summarise(mean_log_dNSAF=mean(log_dNSAF, na.rm=TRUE), .groups="drop") %>% + mutate(col_label=paste(fraction, treatment, sep="_")) %>% + select(MISEVCategory, col_label, mean_log_dNSAF) %>% + pivot_wider(names_from=col_label, values_from=mean_log_dNSAF) %>% + arrange(MISEVCategory) %>% + column_to_rownames("MISEVCategory") %>% + as.matrix() + + svglite(file.path(output_dir, sprintf("marker_category_heatmap_%s.svg", org)), width=8, height=4) + pheatmap(agg_mat, cluster_rows=FALSE, cluster_cols=FALSE, colour=colorRampPalette(c("#88CCEE", "white", "#CC6677"))(50), main=sprintf("Mean log(dNSAF) by MISEV category — %s", org), fontsize=10, border_colour=NA) + dev.off() + + all_marker_tables[[org]] <- marker_table + status <<- record_ok(status, org) + }, error = function(e) { + while (dev.cur() != 1) dev.off() + message(sprintf("EV markers %s: error — %s", org, conditionMessage(e))) + status <<- record_fail(status, org, conditionMessage(e)) + }) } + +write_status(status, output_dir) + # Export .xlsx — one sheet per organism wb <- wb_workbook() for (org in names(all_marker_tables)) { diff --git a/src/comms/r/sections/concordance.R b/src/comms/r/sections/concordance.R index 40dad61..828d61e 100644 --- a/src/comms/r/sections/concordance.R +++ b/src/comms/r/sections/concordance.R @@ -25,6 +25,7 @@ script_dir <- local({ source(file.path(script_dir, "..", "utils", "import.R")) source(file.path(script_dir, "..", "utils", "limma_da.R")) source(file.path(script_dir, "..", "utils", "normalise.R")) +source(file.path(script_dir, "..", "utils", "status.R")) source(file.path(script_dir, "..", "utils", "theme.R")) # Load libraries @@ -51,68 +52,86 @@ lfq_data <- loadLfqFiles(lfq_dir) concordance_stats <- list() organisms <- unique(sample_meta$organism) -concordance_stats <- list() +status <- new_status_tracker("concordance") for (org in organisms) { org_meta <- filter(sample_meta, organism == org) fractions <- unique(org_meta$fraction) for (frac in fractions) { - frac_meta <- filter(org_meta, fraction == frac) - frac_cols <- frac_meta$dnsaf_col - - frac_data <- results_wide %>% - select(proteinId, all_of(frac_cols)) %>% - filter(rowSums(select(., -proteinId) > 0) > 0) - log_mat_dnsaf <- frac_data %>% - column_to_rownames("proteinId") %>% - as.matrix() %>% - logdNSAF() - treatment_vec <- frac_meta %>% - arrange(match(dnsaf_col, frac_cols)) %>% - pull(treatment) - if (length(unique(treatment_vec)) < 2) { - message(sprintf("DA %s %s: fewer than 2 treatment levels in this fraction — skipping", org, frac)) - next - } - da_dnsaf <- runLimmaDA(log_mat_dnsaf, treatment_vec) %>% - classifyDA(lfc_threshold, fdr_threshold) %>% - select(proteinId, log2FC_dNSAF=log2FC, adj_pval_dNSAF=adj_pval, - Abundance_dNSAF=Abundance) - lfq_frac <- filter(lfq_data, Fraction==frac) - if (nrow(lfq_frac) == 0) { - message(sprintf("Concordance %s %s: no LFQ data found — skipping", org, frac)); next - } - - lfq_sample_cols <- setdiff(colnames(lfq_frac), c("proteinId", "Fraction")) - lfq_treatment_vec <- tibble(sample_id = lfq_sample_cols) %>% - inner_join(samples, by = "sample_id") %>% - pull(treatment) - - log_mat_lfq <- lfq_frac %>% - select(proteinId, all_of(lfq_sample_cols)) %>% - column_to_rownames("proteinId") %>% - as.matrix() %>% - log2() - da_lfq <- runLimmaDA(log_mat_lfq, lfq_treatment_vec) %>% - classifyDA(lfc_threshold, fdr_threshold) %>% - select(proteinId, log2FC_LFQ=log2FC, adj_pval_LFQ=adj_pval, Abundance_LFQ=Abundance) - combined <- inner_join(da_dnsaf, da_lfq, by="proteinId") - key <- paste(org, frac, sep = "_") - concordance_stats[[key]] <- combined - r_val <- cor(combined$log2FC_dNSAF, combined$log2FC_LFQ, use = "complete.obs") - scatter <- ggplot(combined, aes(x=log2FC_dNSAF, y=log2FC_LFQ)) + - geom_point(aes(colour=Abundance_dNSAF), alpha=0.6) + - geom_smooth(method="lm", se=FALSE, colour="black", linewidth=0.5) + - scale_colour_manual(values=c("Increased"="#CC6677","Decreased"="#88CCEE","Unchanged"="grey70")) + - theme_comms() + - labs(title=sprintf("LFQ vs dNSAF concordance — %s %s", org, frac), x=expression(log[2](FC)~dNSAF), y=expression(log[2](FC)~LFQ), colour="DA (dNSAF)") + - annotate("text", x=Inf, y=-Inf, hjust=1.1, vjust=-0.5, size=3.5, label=sprintf("r = %.2f (n=%d proteins)", r_val, nrow(combined))) - svglite(file.path(output_dir, sprintf("lfq_vs_dnsaf_%s_%s.svg", frac, org)), width=8, height=7) - print(scatter); dev.off() + tryCatch({ + frac_meta <- filter(org_meta, fraction == frac) + frac_cols <- frac_meta$dnsaf_col + + frac_data <- results_wide %>% + select(proteinId, all_of(frac_cols)) %>% + filter(rowSums(select(., -proteinId) > 0) > 0) + + log_mat_dnsaf <- frac_data %>% + column_to_rownames("proteinId") %>% + as.matrix() %>% + logdNSAF() + treatment_vec <- frac_meta %>% + arrange(match(dnsaf_col, frac_cols)) %>% + pull(treatment) + if (length(unique(treatment_vec)) < 2) { + reason <- sprintf("fewer than 2 treatment levels in fraction %s", frac) + message(sprintf("Concordance %s %s: %s — skipping", org, frac, reason)) + status <<- record_skip(status, org, reason) + next + } + da_dnsaf <- runLimmaDA(log_mat_dnsaf, treatment_vec) %>% + classifyDA(lfc_threshold, fdr_threshold) %>% + select(proteinId, log2FC_dNSAF=log2FC, adj_pval_dNSAF=adj_pval, Abundance_dNSAF=Abundance) + + lfq_frac <- filter(lfq_data, Fraction == frac) + if (nrow(lfq_frac) == 0) { + reason <- sprintf("no LFQ data for fraction %s", frac) + message(sprintf("Concordance %s %s: %s — skipping", org, frac, reason)) + status <<- record_skip(status, org, reason) + next + } + + lfq_sample_cols <- setdiff(colnames(lfq_frac), c("proteinId", "Fraction")) + lfq_treatment_vec <- tibble(sample_id=lfq_sample_cols) %>% + inner_join(samples, by="sample_id") %>% + pull(treatment) + + log_mat_lfq <- lfq_frac %>% + select(proteinId, all_of(lfq_sample_cols)) %>% + column_to_rownames("proteinId") %>% + as.matrix() %>% + log2() + da_lfq <- runLimmaDA(log_mat_lfq, lfq_treatment_vec) %>% + classifyDA(lfc_threshold, fdr_threshold) %>% + select(proteinId, log2FC_LFQ=log2FC, adj_pval_LFQ=adj_pval, Abundance_LFQ=Abundance) + + combined <- inner_join(da_dnsaf, da_lfq, by="proteinId") + key <- paste(org, frac, sep = "_") + concordance_stats[[key]] <- combined + + r_val <- cor(combined$log2FC_dNSAF, combined$log2FC_LFQ, use="complete.obs") + scatter <- ggplot(combined, aes(x=log2FC_dNSAF, y=log2FC_LFQ)) + + geom_point(aes(colour=Abundance_dNSAF), alpha=0.6) + + geom_smooth(method="lm", se=FALSE, colour="black", linewidth=0.5) + + scale_colour_manual(values=c("Increased"="#CC6677", "Decreased"="#88CCEE", "Unchanged"="grey70")) + + theme_comms() + + labs(title=sprintf("LFQ vs dNSAF concordance — %s %s", org, frac), x=expression(log[2](FC)~dNSAF), y=expression(log[2](FC)~LFQ), colour="DA (dNSAF)") + + annotate("text", x=Inf, y =-Inf, hjust=1.1, vjust=-0.5, size=3.5, label=sprintf("r = %.2f (n=%d proteins)", r_val, nrow(combined))) + svglite(file.path(output_dir, sprintf("lfq_vs_dnsaf_%s_%s.svg", frac, org)), width=8, height=7) + print(scatter); dev.off() + + status <<- record_ok(status, org) + }, error = function(e) { + while (dev.cur() != 1) dev.off() + message(sprintf("Concordance %s %s: error — %s", org, frac, conditionMessage(e))) + status <<- record_fail(status, org, conditionMessage(e)) + }) } } +write_status(status, output_dir) + # Export .xlsx spreadsheet wb <- wb_workbook() for (key in names(concordance_stats)) { diff --git a/src/comms/r/sections/da.R b/src/comms/r/sections/da.R index 8adc634..6dffbdc 100644 --- a/src/comms/r/sections/da.R +++ b/src/comms/r/sections/da.R @@ -12,6 +12,7 @@ organism_prefix <- args[6] min_reps <- as.integer(args[7]) lfc_threshold <- as.numeric(args[8]) fdr_threshold <- as.numeric(args[9]) +top_n <- as.integer(args[10]) # Get script directory for path traversal script_dir <- local({ @@ -24,6 +25,7 @@ script_dir <- local({ source(file.path(script_dir, "..", "utils", "import.R")) source(file.path(script_dir, "..", "utils", "limma_da.R")) source(file.path(script_dir, "..", "utils", "normalise.R")) +source(file.path(script_dir, "..", "utils", "status.R")) source(file.path(script_dir, "..", "utils", "theme.R")) # Load libraries @@ -47,76 +49,91 @@ if (length(treatments) != 2) stop("DA section requires exactly two treatment lev # Initialise list for differentially abundant results da_results_all <- list() +status <- new_status_tracker("da") + for (org in organisms) { org_meta <- filter(sample_meta, organism == org) fractions <- unique(org_meta$fraction) for (frac in fractions) { - frac_meta <- filter(org_meta, fraction==frac) - frac_cols <- frac_meta$dnsaf_col - - frac_data <- results_wide %>% - select(proteinId, proteinAnnotation, all_of(frac_cols)) %>% - filter(rowSums(select(., all_of(frac_cols)) > 0) > 0) - - for (trt in treatments) { - trt_cols <- filter(frac_meta, treatment==trt)$dnsaf_col - frac_data[[paste0("n_", trt)]] <- rowSums(select(frac_data, all_of(trt_cols)) > 0) - } - frac_data <- filter(frac_data, if_any(starts_with("n_"), ~. >= min_reps)) - - if (nrow(frac_data) < 5) { - message(sprintf("DA %s %s: too few proteins after replicate filter (%d) — skipping", org, frac, nrow(frac_data))); next - } - - log_mat <- frac_data %>% - select(proteinId, all_of(frac_cols)) %>% - column_to_rownames("proteinId") %>% - as.matrix() %>% - logdNSAF() - treatment_vec <- frac_meta %>% - arrange(match(dnsaf_col, frac_cols)) %>% - pull(treatment) - - if (length(unique(treatment_vec)) < 2) { - message(sprintf("DA %s %s: fewer than 2 treatment levels in this fraction — skipping", org, frac)) - next - } - - da_res <- runLimmaDA(log_mat, treatment_vec) %>% - classifyDA(lfc_threshold, fdr_threshold) %>% - left_join(select(frac_data, proteinId, proteinAnnotation), by="proteinId") - - key <- paste(org, frac, sep="_") - da_results_all[[key]] <- da_res - - # Generate volcano plot - top_labels <- filter(da_res, Abundance != "Unchanged") %>% slice_min(adj_pval, n=20) - volcano <- ggplot(da_res, aes(x=log2FC, y=-log10(adj_pval), colour=Abundance)) + - geom_point(alpha=0.7, size=1.5) + - geom_hline(yintercept=-log10(fdr_threshold), linetype="dashed", colour="grey50") + - geom_vline(xintercept=c(-lfc_threshold, lfc_threshold), linetype="dashed", colour="grey50") + - geom_text_repel(data=top_labels, aes(label=proteinAnnotation), size=3, max.overlaps=15) + - scale_colour_manual(values=c("Increased"="#CC6677","Decreased"="#88CCEE","Unchanged"="grey70")) + - theme_comms() + - labs(title=sprintf("DA — %s %s (%s vs %s)", org, frac, treatments[2], treatments[1]), x=expression(log[2](FC)), y=expression(-log[10](adj.p))) - svglite(file.path(output_dir, sprintf("volcano_%s_%s.svg", frac, org)), width=10, height=7) - print(volcano); dev.off() + tryCatch({ + frac_meta <- filter(org_meta, fraction == frac) + frac_cols <- frac_meta$dnsaf_col + + frac_data <- results_wide %>% + select(proteinId, proteinAnnotation, all_of(frac_cols)) %>% + filter(rowSums(select(., all_of(frac_cols)) > 0) > 0) + + for (trt in treatments) { + trt_cols <- filter(frac_meta, treatment == trt)$dnsaf_col + frac_data[[paste0("n_", trt)]] <- rowSums(select(frac_data, all_of(trt_cols)) > 0) + } + frac_data <- filter(frac_data, if_any(starts_with("n_"), ~. >= min_reps)) + + if (nrow(frac_data) < 5) { + reason <- sprintf("too few proteins after replicate filter (%d) in fraction %s", nrow(frac_data), frac) + message(sprintf("DA %s %s: %s — skipping", org, frac, reason)) + status <<- record_skip(status, org, reason) + next + } + + log_mat <- frac_data %>% + select(proteinId, all_of(frac_cols)) %>% + column_to_rownames("proteinId") %>% + as.matrix() %>% + logdNSAF() + treatment_vec <- frac_meta %>% + arrange(match(dnsaf_col, frac_cols)) %>% + pull(treatment) + + if (length(unique(treatment_vec)) < 2) { + reason <- sprintf("fewer than 2 treatment levels in fraction %s", frac) + message(sprintf("DA %s %s: %s — skipping", org, frac, reason)) + status <<- record_skip(status, org, reason) + next + } + + da_res <- runLimmaDA(log_mat, treatment_vec) %>% + classifyDA(lfc_threshold, fdr_threshold) %>% + left_join(select(frac_data, proteinId, proteinAnnotation), by = "proteinId") + + key <- paste(org, frac, sep = "_") + da_results_all[[key]] <- da_res + + top_labels <- filter(da_res, Abundance != "Unchanged") %>% + slice_min(adj_pval, n=top_n) + volcano <- ggplot(da_res, aes(x=log2FC, y=-log10(adj_pval), colour=Abundance)) + + geom_point(alpha=0.7, size=1.5) + + geom_hline(yintercept=-log10(fdr_threshold), linetype="dashed", colour="grey50") + + geom_vline(xintercept=c(-lfc_threshold, lfc_threshold), linetype="dashed", colour="grey50") + + geom_text_repel(data=top_labels, aes(label=proteinAnnotation), size=3, max.overlaps=15) + + scale_colour_manual(values=c("Increased"="#CC6677", "Decreased"="#88CCEE", "Unchanged"="grey70")) + + theme_comms() + + labs(title=sprintf("DA — %s %s (%s vs %s)", org, frac, treatments[2], treatments[1]), x=expression(log[2](FC)), y=expression(-log[10](adj.p))) + svglite(file.path(output_dir, sprintf("volcano_%s_%s.svg", frac, org)), width=10, height=7) + print(volcano); dev.off() + + status <<- record_ok(status, org) + }, error = function(e) { + while (dev.cur() != 1) dev.off() + message(sprintf("DA %s %s: error — %s", org, frac, conditionMessage(e))) + status <<- record_fail(status, org, conditionMessage(e)) + }) } } +write_status(status, output_dir) + # Venn diagrams per organism for (org in organisms) { org_results <- da_results_all[str_starts(names(da_results_all), org)] - da_up_sets <- lapply(org_results, function(x) filter(x, Abundance == "Increased")$proteinId) + da_up_sets <- lapply(org_results, function(x) filter(x, Abundance == "Increased")$proteinId) da_down_sets <- lapply(org_results, function(x) filter(x, Abundance == "Decreased")$proteinId) names(da_up_sets) <- names(da_down_sets) <- str_remove(names(org_results), paste0(org, "_")) if (length(da_up_sets) >= 2) { - venn_up <- venn.diagram(da_up_sets, filename=NULL, disable.logging=TRUE, - category.names=names(da_up_sets)) + venn_up <- venn.diagram(da_up_sets, filename=NULL, disable.logging=TRUE, category.names=names(da_up_sets)) ggsave(file.path(output_dir, sprintf("venn_da_up_%s.svg", org)), venn_up) - venn_down <- venn.diagram(da_down_sets, filename=NULL, disable.logging=TRUE, - category.names=names(da_down_sets)) + venn_down <- venn.diagram(da_down_sets, filename=NULL, disable.logging=TRUE, category.names=names(da_down_sets)) ggsave(file.path(output_dir, sprintf("venn_da_down_%s.svg", org)), venn_down) } } diff --git a/src/comms/r/sections/pca.R b/src/comms/r/sections/pca.R index 4517c75..5797811 100644 --- a/src/comms/r/sections/pca.R +++ b/src/comms/r/sections/pca.R @@ -20,6 +20,7 @@ script_dir <- local({ # Import utility functions source(file.path(script_dir, "..", "utils", "import.R")) +source(file.path(script_dir, "..", "utils", "status.R")) source(file.path(script_dir, "..", "utils", "theme.R")) # Load libraries @@ -28,6 +29,9 @@ library(ggfortify) library(ggrepel) library(svglite) +# Define variable for minimum samples required to run clustering +MIN_SAMPLES_FOR_PCA <- max(min_reps, 2) + # Import files ref_info <- loadRefInfo(ref_info_path) cont_info <- loadContInfo(cont_csv_path) @@ -38,36 +42,55 @@ dnsaf_cols <- colnames(results_wide)[startsWith(colnames(results_wide), "dNSAF_" sample_meta <- buildSampleMetadata(str_remove(dnsaf_cols, "dNSAF_"), samples) organisms <- unique(sample_meta$organism) +status <- new_status_tracker("pca") + for (org in organisms) { - org_meta <- filter(sample_meta, organism == org) - org_cols <- org_meta$dnsaf_col - label_map <- setNames(org_meta$sample_id, org_meta$dnsaf_col) + tryCatch({ + org_meta <- filter(sample_meta, organism == org) + org_cols <- org_meta$dnsaf_col + label_map <- setNames(org_meta$sample_id, org_meta$dnsaf_col) + + if (length(org_cols) < MIN_SAMPLES_FOR_PCA) { + reason <- sprintf("only %d sample(s) available (need >= %d for clustering)", length(org_cols), MIN_SAMPLES_FOR_PCA) + message(sprintf("PCA %s: %s - skipping", org, reason)) + status <<- record_skip(status, org, reason) + next + } - org_data <- results_wide %>% - filter(rowSums(select(., all_of(org_cols)) > 0) > 0) + org_data <- results_wide %>% + filter(rowSums(select(., all_of(org_cols)) > 0) > 0) - pca_mat <- org_data %>% - select(proteinId, all_of(org_cols)) %>% - column_to_rownames("proteinId") %>% - rename_with(~label_map[.]) %>% - t() + pca_mat <- org_data %>% + select(proteinId, all_of(org_cols)) %>% + column_to_rownames("proteinId") %>% + rename_with(~label_map[.]) %>% + t() - k <- max(2, min(length(unique(org_meta$fraction)), nrow(pca_mat) - 1)) - pca_data <- clara(pca_mat, k=k, metric="euclidean", stand=FALSE, samples=500, sampsize=nrow(pca_mat), pamLike=TRUE, correct.d=TRUE) + k <- max(2, min(length(unique(org_meta$fraction)), nrow(pca_mat) - 1)) + pca_data <- clara(pca_mat, k=k, metric="euclidean", stand=FALSE, samples=500, sampsize=nrow(pca_mat), pamLike=TRUE, correct.d=TRUE) - pca_plot <- autoplot(pca_data, frame=TRUE, frame.type="t", size=5) + - theme_comms() + - geom_text_repel(label=rownames(pca_mat), size=4, box.padding=0.5, point.padding=0.75, direction="both", force=15, max.overlaps=Inf) + - scale_color_manual(values=COMMS_COLOURS) + - scale_fill_manual(values=COMMS_COLOURS) + - labs(colour="Cluster", fill="Cluster", title=sprintf("PCA — %s", org)) - svglite(file.path(output_dir, sprintf("pca_%s.svg", org)), width=10, height=8) - print(pca_plot); dev.off() + pca_plot <- autoplot(pca_data, frame=TRUE, frame.type="t", size=5) + + theme_comms() + + geom_text_repel(label=rownames(pca_mat), size=4, box.padding=0.5, point.padding=0.75, direction="both", force=15, max.overlaps=Inf) + + scale_color_manual(values = COMMS_COLOURS) + + scale_fill_manual(values = COMMS_COLOURS) + + labs(colour="Cluster", fill="Cluster", title=sprintf("PCA — %s", org)) + svglite(file.path(output_dir, sprintf("pca_%s.svg", org)), width=10, height=8) + print(pca_plot); dev.off() - dist_mat <- dist(scale(pca_mat), method="euclidean") - hc <- hclust(dist_mat, method="average") - svglite(file.path(output_dir, sprintf("dendrogram_%s.svg", org)), width=10, height=6) - plot(hc, main=sprintf("Sample clustering — %s", org), xlab="", sub="", ylab="Distance", cex=0.9) - dev.off() + dist_mat <- dist(scale(pca_mat), method="euclidean") + hc <- hclust(dist_mat, method="average") + svglite(file.path(output_dir, sprintf("dendrogram_%s.svg", org)), width=10, height=6) + plot(hc, main=sprintf("Sample clustering — %s", org), xlab="", sub="", ylab="Distance", cex=0.9) + dev.off() + + status <<- record_ok(status, org) + }, error = function(e) { + while (dev.cur() != 1) dev.off() + message(sprintf("PCA %s: error — %s", org, conditionMessage(e))) + status <<- record_fail(status, org, conditionMessage(e)) + }) } + +write_status(status, output_dir) message("PCA section complete") \ No newline at end of file diff --git a/src/comms/r/sections/qc.R b/src/comms/r/sections/qc.R index 2c303ba..ab625d9 100644 --- a/src/comms/r/sections/qc.R +++ b/src/comms/r/sections/qc.R @@ -21,6 +21,7 @@ script_dir <- local({ # Import utility functions source(file.path(script_dir, "..", "utils", "import.R")) source(file.path(script_dir, "..", "utils", "normalise.R")) +source(file.path(script_dir, "..", "utils", "status.R")) source(file.path(script_dir, "..", "utils", "theme.R")) # Load libraries @@ -38,63 +39,89 @@ dnsaf_cols <- colnames(results_wide)[startsWith(colnames(results_wide), "dNSAF_" sample_meta <- buildSampleMetadata(str_remove(dnsaf_cols, "dNSAF_"), samples) organisms <- unique(sample_meta$organism) +status <- new_status_tracker("qc") + for (org in organisms) { - org_meta <- filter(sample_meta, organism == org) - org_cols <- org_meta$dnsaf_col - label_map <- setNames(org_meta$sample_id, org_meta$dnsaf_col) - - org_data <- results_wide %>% - filter(rowSums(select(., all_of(org_cols)) > 0) > 0) - - # Per-sample dNSAF density plot - dnsaf_long <- org_data %>% - select(proteinId, all_of(org_cols)) %>% - pivot_longer(-proteinId, names_to="Sample", values_to="dNSAF") %>% - filter(dNSAF > 0) %>% - mutate(log_dNSAF=log(dNSAF), Sample=label_map[Sample]) - density_plot <- ggplot(dnsaf_long, aes(x=log_dNSAF, colour=Sample)) + - geom_density() + theme_comms() + - labs(x="log(dNSAF)", y="Density", title=sprintf("Per-sample dNSAF distributions — %s", org)) + - theme(legend.position="bottom") - svglite(file.path(output_dir, sprintf("dnsaf_distributions_%s.svg", org)), width=10, height=6) - print(density_plot); dev.off() - - # Total spectral counts per sample - spec_counts <- bind_rows(lapply(org_cols, function(col) { - nm <- str_remove(col, "dNSAF_") - if (!nm %in% names(results_list)) return(NULL) - tibble(Sample=label_map[col], TotalSpectra=sum(results_list[[nm]]$`RAW`, na.rm=TRUE)) - })) %>% compact() %>% bind_rows() - counts_plot <- ggplot(spec_counts, aes(x=Sample, y=TotalSpectra)) + - geom_col(fill="#88CCEE") + theme_comms() + - theme(axis.text.x=element_text(angle=45, hjust=1)) + - labs(x=NULL, y="Total spectral counts", title=sprintf("Spectral counts per sample — %s", org)) - svglite(file.path(output_dir, sprintf("spectral_counts_per_sample_%s.svg", org)), width=10, height=5) - print(counts_plot); dev.off() - - # Missing-value upset plot - presence_matrix <- org_data %>% - select(all_of(org_cols)) %>% - mutate(across(everything(), ~as.integer(. > 0))) - colnames(presence_matrix) <- label_map[colnames(presence_matrix)] - svglite(file.path(output_dir, sprintf("missing_values_upset_%s.svg", org)), width=12, height=7) - upset(as.data.frame(presence_matrix), nsets=ncol(presence_matrix), order.by="freq", mainbar.y.label="Proteins", sets.x.label="Proteins detected") - dev.off() - - # Presence/absence heatmap - svglite(file.path(output_dir, sprintf("presence_absence_heatmap_%s.svg", org)), width=10, height=8) - pheatmap(as.matrix(presence_matrix), color=c("white", "#117733"), legend_breaks=c(0, 1), legend_labels=c("Absent", "Present"), main=sprintf("Protein presence/absence — %s", org), fontsize=10) - dev.off() - - # QC summary Excel - n_detected <- org_data %>% - summarise(across(all_of(org_cols), ~sum(. > 0))) %>% - pivot_longer(everything(), names_to="dnsaf_col", values_to="ProteinsDetected") %>% - mutate(Sample=label_map[dnsaf_col]) %>% - select(Sample, ProteinsDetected) - qc_summary <- left_join(spec_counts, n_detected, by="Sample") - wb <- wb_workbook() - wb$add_worksheet(org); wb$add_data(org, qc_summary) - wb_save(wb, file.path(output_dir, sprintf("qc_summary_%s.xlsx", org))) + tryCatch({ + org_meta <- filter(sample_meta, organism == org) + org_cols <- org_meta$dnsaf_col + label_map <- setNames(org_meta$sample_id, org_meta$dnsaf_col) + + org_data <- results_wide %>% + filter(rowSums(select(., all_of(org_cols)) > 0) > 0) + + if (nrow(org_data) == 0) { + reason <- "no proteins detected for this organism" + message(sprintf("QC %s: %s — skipping", org, reason)) + status <<- record_skip(status, org, reason) + next + } + + # Per-sample dNSAF density plot (works fine with a single sample) + dnsaf_long <- org_data %>% + select(proteinId, all_of(org_cols)) %>% + pivot_longer(-proteinId, names_to="Sample", values_to="dNSAF") %>% + filter(dNSAF > 0) %>% + mutate(log_dNSAF=log(dNSAF), Sample=label_map[Sample]) + density_plot <- ggplot(dnsaf_long, aes(x=log_dNSAF, colour=Sample)) + + geom_density() + theme_comms() + + labs(x="log(dNSAF)", y="Density", title=sprintf("Per-sample dNSAF distributions — %s", org)) + + theme(legend.position="bottom") + svglite(file.path(output_dir, sprintf("dnsaf_distributions_%s.svg", org)), width=10, height=6) + print(density_plot); dev.off() + + # Total spectral counts per sample (also fine with a single sample) + spec_counts <- bind_rows(lapply(org_cols, function(col) { + nm <- str_remove(col, "dNSAF_") + if (!nm %in% names(results_list)) return(NULL) + tibble(Sample=label_map[col], TotalSpectra=sum(results_list[[nm]]$`RAW`, na.rm=TRUE)) + })) %>% compact() %>% bind_rows() + counts_plot <- ggplot(spec_counts, aes(x=Sample, y=TotalSpectra)) + + geom_col(fill="#88CCEE") + theme_comms() + + theme(axis.text.x=element_text(angle=45, hjust=1)) + + labs(x=NULL, y="Total spectral counts", title=sprintf("Spectral counts per sample — %s", org)) + svglite(file.path(output_dir, sprintf("spectral_counts_per_sample_%s.svg", org)), width=10, height=5) + print(counts_plot); dev.off() + + # Missing-value upset plot and presence/absence heatmap need >= 2 samples to mean anything + if (length(org_cols) >= 2) { + presence_matrix <- org_data %>% + select(all_of(org_cols)) %>% + mutate(across(everything(), ~as.integer(. > 0))) + colnames(presence_matrix) <- label_map[colnames(presence_matrix)] + + svglite(file.path(output_dir, sprintf("missing_values_upset_%s.svg", org)), width = 12, height = 7) + upset(as.data.frame(presence_matrix), nsets = ncol(presence_matrix), order.by = "freq", + mainbar.y.label = "Proteins", sets.x.label = "Proteins detected") + dev.off() + + svglite(file.path(output_dir, sprintf("presence_absence_heatmap_%s.svg", org)), width = 10, height = 8) + pheatmap(as.matrix(presence_matrix), color = c("white", "#117733"), legend_breaks = c(0, 1), + legend_labels = c("Absent", "Present"), + main = sprintf("Protein presence/absence — %s", org), fontsize = 10) + dev.off() + } else { + message(sprintf("QC %s: only 1 sample — skipping upset plot and presence/absence heatmap", org)) + } + + # QC summary Excel + n_detected <- org_data %>% + summarise(across(all_of(org_cols), ~sum(. > 0))) %>% + pivot_longer(everything(), names_to = "dnsaf_col", values_to = "ProteinsDetected") %>% + mutate(Sample = label_map[dnsaf_col]) %>% + select(Sample, ProteinsDetected) + qc_summary <- left_join(spec_counts, n_detected, by = "Sample") + wb <- wb_workbook() + wb$add_worksheet(org); wb$add_data(org, qc_summary) + wb_save(wb, file.path(output_dir, sprintf("qc_summary_%s.xlsx", org))) + + status <<- record_ok(status, org) + }, error = function(e) { + while (dev.cur() != 1) dev.off() + message(sprintf("QC %s: error — %s", org, conditionMessage(e))) + status <<- record_fail(status, org, conditionMessage(e)) + }) } + +write_status(status, output_dir) message("QC section complete") \ No newline at end of file diff --git a/src/comms/r/utils/status.R b/src/comms/r/utils/status.R new file mode 100644 index 0000000..ebe31da --- /dev/null +++ b/src/comms/r/utils/status.R @@ -0,0 +1,48 @@ +#!/bin/R +# status.R: per-organism status tracking utilities + +library(jsonlite) + +# new_status_tracker: initialise an empty tracker for a section +new_status_tracker <- function(section) { + list(section = section, organisms = list(), reasons = list()) +} + +# record_ok: mark organism as having produced at least one successful output +record_ok <- function(tracker, organism) { + current <- tracker$organisms[[organism]] + if (is.null(current) || current == "skipped") { + tracker$organisms[[organism]] <- "ok" + } + tracker +} + +# record_skip: mark organism as skipped for this attempt (insufficient data) +record_skip <- function(tracker, organism, reason) { + if (is.null(tracker$organisms[[organism]])) { + tracker$organisms[[organism]] <- "skipped" + } + tracker$reasons[[organism]] <- c(tracker$reasons[[organism]], reason) + tracker +} + +# record_fail: mark organism as failed (unhandled error) +record_fail <- function(tracker, organism, reason) { + tracker$organisms[[organism]] <- "failed" + tracker$reasons[[organism]] <- c(tracker$reasons[[organism]], reason) + tracker +} + +# write_status: serialise a tracker to /_status.json +write_status <- function(tracker, output_dir) { + reasons_flat <- lapply(tracker$reasons, function(r) paste(r, collapse = "; ")) + payload <- list( + section = tracker$section, + organisms = tracker$organisms, + reasons = if (length(reasons_flat)) reasons_flat else setNames(list(), character(0)) + ) + writeLines( + toJSON(payload, auto_unbox = TRUE, pretty = TRUE), + file.path(output_dir, "_status.json") + ) +} \ No newline at end of file diff --git a/src/comms/utils/crux.py b/src/comms/utils/crux.py index 53ffcb2..9aa7886 100644 --- a/src/comms/utils/crux.py +++ b/src/comms/utils/crux.py @@ -114,20 +114,17 @@ def paramMedic(crux_bin: Path, mzml_file: Path, out_dir: Path) -> bool: return runCrux(crux_bin, 'param-medic', args) # -- tideSearch: returns True if Tide-search completed successfully for the given mzML file, False on failure -def tideSearch(crux_bin: Path, mzml_file: Path, index_dir: Path, out_dir: Path, fileroot: str, config: dict, threads, precursor_tol=None, mz_bin_width=None) -> bool: +def tideSearch(crux_bin: Path, mzml_file: Path, index_dir: Path, out_dir: Path, fileroot: str, config: dict) -> bool: logMsg.debug(f'tide-search: {mzml_file.name}') - prec = precursor_tol or config['search']['precursor_tolerance_ppm'] - bin_width = mz_bin_width or config['search']['mz_bin_width'] - logMsg.debug(f'Precursor tolerance {prec} ppm, m/z bin width {bin_width}') args = [ '--verbosity', '40', - '--num-threads', threads, + '--num-threads', config['search']['threads'], '--spectrum-parser', 'pwiz', - '--precursor-window', str(prec), + '--precursor-window', str(config['search']['precursor_tolerance_ppm']), '--precursor-window-type', 'ppm', - '--mz-bin-width', str(bin_width), + '--mz-bin-width', str(config['search']['mz_bin_width']), '--score-function', config['search']['score_function'], - '--min-peaks', str(config['search']['min_peaks']), + '--min-peaks', str(config['search']['min_peaks'],), '--missed-cleavages', str(config['index']['missed_cleavages']), '--output-dir', str(out_dir), '--fileroot', fileroot, @@ -142,13 +139,13 @@ def percolator(crux_bin: Path, target_psm_file: Path, database: Path, out_dir: P logMsg.debug(f'percolator: {target_psm_file.name}') args = [ '--verbosity', '40', - '--protein-enzyme', config['percolator']['protein_enzyme'], + '--protein-enzyme', config['rescore']['protein_enzyme'], '--output-dir', str(out_dir), '--fileroot', fileroot, '--overwrite', 'T', str(target_psm_file), ] - if config['percolator']['picked_protein']: + if config['rescore']['picked_protein']: args = ['--picked-protein', str(database)] + args return runCrux(crux_bin, 'percolator', args) @@ -230,8 +227,7 @@ def lfq(crux_bin, psm_files, mzml_files, out_dir, fileroot, config) -> bool: return ok # -- _tomlToCrux: helper function returning 'T' if True and 'F' if False -def _tomlToCrux(val: bool): - if val: - return 'T' - else: - return 'F' \ No newline at end of file +def _tomlToCrux(val) -> str: + if isinstance(val, str): + val = val.strip().lower() == 'true' + return 'T' if val else 'F' \ No newline at end of file diff --git a/src/comms/utils/installrdeps.py b/src/comms/utils/installrdeps.py new file mode 100644 index 0000000..fdbc50c --- /dev/null +++ b/src/comms/utils/installrdeps.py @@ -0,0 +1,96 @@ +''' +comMS R dependency wrapper utilities +''' + +# -- Import external dependencies +import json, shutil, subprocess +from importlib.resources import files as pkg_files +from pathlib import Path +from rich import print + +# -- Import internal functions +from comms.utils.log import logMsg + +# -- _r_script: returns Path to a script under comms/r/ +def _r_script(name: str) -> Path: + return pkg_files('comms').joinpath(f'r/{name}') + +# -- _print_dependency_table: returns None but prints a table to output which lists installed/unavailable dependencies +def _print_dependency_table(deps_dict: dict[str, list[str]]) -> None: + installed_deps = deps_dict['installed'] + missing_deps = deps_dict['missing'] + if len(installed_deps) > 0: + logMsg.info(f'[bold green]✓ Installed packages ({len(installed_deps)})[/bold green]: {", ".join(installed_deps)}') + if len(missing_deps) > 0: + logMsg.info(f'[bold red]✗ Missing packages ({len(missing_deps)})[/bold red]: {", ".join(missing_deps)}') + logMsg.info('To install missing packages, run [bold]comms r-utils install[/bold]') + +# -- check_r_dependencies: returns {'installed': [...], 'missing': [...]}, or None if Rscript itself isn't callable / the check failed +def check_r_dependencies(rscript: str = 'Rscript') -> dict[str, list[str]] | None: + if shutil.which(rscript) is None: + logMsg.error(f'Rscript not callable: {rscript}') + return None + logMsg.debug(f'R available at {rscript}') + script_path = _r_script('/deps/check_deps.R') + result = subprocess.run([rscript, '--vanilla', str(script_path)], capture_output=True, text=True) + logMsg.debug(f'R dependency check ran successfully') + if result.returncode != 0: + logMsg.warn(f'Dependency check failed: {result.stderr}') + return None + # Parse returned result as JSON + try: + parsed = json.loads(result.stdout.strip()) + _print_dependency_table(parsed) + return parsed + except Exception as e: + logMsg.error(f'Could not parse depencency check output: {e}') + return None + +# -- install_r_dependencies: runs install_deps.R, streaming its messages through logMsg.info; returns True on success (including "nothing to do") +def install_r_dependencies(rscript: str = 'Rscript') -> bool: + if shutil.which(rscript) is None: + logMsg.error(f'Rscript not callable: {rscript}') + return False + script_path = _r_script('deps/install_deps.R') + logMsg.info('Installing R report dependencies (this can take a few minutes)...') + process = subprocess.Popen( + [rscript, '--vanilla', str(script_path)], + stdout=subprocess.PIPE, stderr=subprocess.STDOUT, text=True, + ) + for line in process.stdout: + line = line.rstrip() + if line: + logMsg.info(line) + process.wait() + if process.returncode != 0: + logMsg.error('R dependency installation failed; see above output') + return False + logMsg.info('All R report dependencies installed') + return True + +# -- install_r_dependencies_terminal: wrapper for install_r_dependencies to print confirm installation +def install_r_dependencies_terminal(rscript: str = 'Rscript') -> None: + # Check that Rscript is available + if shutil.which(rscript) is None: + logMsg.error(f'Rscript not callable: {rscript}') + raise SystemExit(1) + # Check if anything needs installing + try: + script_path = _r_script('/deps/check_deps.R') + result = subprocess.run([rscript, '--vanilla', str(script_path)], capture_output=True, text=True) + parsed = json.loads(result.stdout.strip()) + missing = parsed['missing'] + if len(missing) > 0: + logMsg.info(f'{len(missing)} {"dependencies need" if len(missing) > 1 else "dependency needs"} to be installed: {", ".join(d for d in missing)}') + while True: + user_confirmation = logMsg.input('Install these packages? [dim](y/N)[/dim]', choices=['y','n'], default='n', case_sensitive=False, show_default=False, show_choices=False).lower() + if user_confirmation in ['', 'n']: + logMsg.info(f'Cancelled dependency installation') + break + if user_confirmation == 'y': + install_r_dependencies() + break + else: + logMsg.info(f'No dependencies missing') + except Exception as e: + logMsg.error(f'Error while checking for missing dependencies: {e}') \ No newline at end of file diff --git a/src/comms/utils/log.py b/src/comms/utils/log.py index 5e4d48f..8c80bd5 100644 --- a/src/comms/utils/log.py +++ b/src/comms/utils/log.py @@ -2,10 +2,11 @@ Shared utility functions: logging ''' # -- Import external dependencies -import atexit, logging, shutil, sys, tempfile +import atexit, logging, shutil, sys, tempfile, time from pathlib import Path from rich.console import Console, ConsoleRenderable from rich.logging import RichHandler +from rich.prompt import Prompt from rich.text import Text # -- Define level colours for RichHandler @@ -13,6 +14,7 @@ 'DEBUG': 'color(67)', 'PROGRESS': 'color(75)', 'INFO': 'color(33)', + 'INPUT': 'color(28)', 'WARNING': 'color(178)', 'ERROR': 'color(160)', 'CRITICAL': 'color(124)', @@ -28,10 +30,23 @@ def _progress(self, message, *args, **kwargs): logging.Logger.progress = _progress +# -- Register custom INPUT logging level +INPUT = 45 +logging.addLevelName(INPUT, 'INPUT') + +def _input(self, message, *args, **kwargs): + if self.isEnabledFor(INPUT): + self._log(INPUT, message, args, **kwargs) + +logging.Logger.input = _input + + # -- Define custom RichHandler subclass (CommsRichHandler) to allow custom formatting class CommsRichHandler(RichHandler): def emit(self, record: logging.LogRecord) -> None: log_state._emitted = True + if getattr(record, '_suppress_console', False): + return super().emit(record) def render_message(self, record: logging.LogRecord, message: str) -> 'ConsoleRenderable': level_colour = _LEVEL_COLOURS.get(record.levelname, 'white') @@ -61,6 +76,46 @@ def info(cls, msg: str): if cls._instance: cls._instance.logger.info(msg) @classmethod + def input(cls, msg: str, **prompt_kwargs) -> str | None: + if not cls._instance: + return None + logger = cls._instance.logger + if not logger.isEnabledFor(INPUT): + return None + interactive = sys.stdin.isatty() + if interactive: + level_colour = _LEVEL_COLOURS.get('INPUT', 'white') + console = Console(stderr=True) + timestamp = time.strftime('%Y-%m-%d %H:%M:%S') + styled_prompt = ( + f"[dim]{timestamp}[/dim] | " + f"[bold]{logger.name}[/bold] | " + f"[bold {level_colour}]INPUT[/] | [white]{msg}[/]" + ) + prompt_obj = Prompt( + styled_prompt, + console=console, + choices=prompt_kwargs.get('choices'), + show_default=prompt_kwargs.get('show_default', True), + show_choices=prompt_kwargs.get('show_choices', True), + ) + prompt_obj.case_sensitive = prompt_kwargs.get('case_sensitive', True) + default = prompt_kwargs.get('default', ...) + answer = prompt_obj(default=default, stream=prompt_kwargs.get('stream')) + # Work out how many terminal rows the prompt (and the typed answer, if echoed) actually occupied, so wrapped prompts get fully erased rather than leaving fragments behind + rendered = prompt_obj.make_prompt(default) + total_len = len(rendered.plain) + len(str(answer)) + width = console.width or 80 + rows = max(1, -(-total_len // width)) # ceil division + for _ in range(rows): + console.file.write("\x1b[1A\x1b[2K") + console.file.write("\r") + console.file.flush() + else: + answer = prompt_kwargs.get('default') + logger.input(f"{msg}: {answer}") + return answer + @classmethod def warn(cls, msg: str): if cls._instance: cls._instance.logger.warning(msg) @@ -77,6 +132,7 @@ class LogState: _file_handler: logging.FileHandler | None = None _emitted: bool = False _atexit_registered: bool = False + _pipeline_log_paths: list[Path] | None = None # -- Define custom logging.Formatter subclas (PlainFormatter) to strip Rich markup before writing to file class PlainFormatter(logging.Formatter): @@ -147,6 +203,31 @@ def configureFileLogging(out_dir: Path): logger.addHandler(handler) log_state._file_handler = handler _removeTempLog() + # Track this log file if pipeline is aggregating logs + if log_state._pipeline_log_paths is not None: + log_state._pipeline_log_paths.append(final_path) + +# -- startPipelineLogging: begin tracking per-command log files for later aggregation into pipeline.log +def startPipelineLogging(): + log_state._pipeline_log_paths = [] + +# -- concatenatePipelineLog: concatenate per-command log files into pipeline.log +def concatenatePipelineLog(out_path: Path) -> Path | None: + paths = log_state._pipeline_log_paths or [] + log_state._pipeline_log_paths = None # stop tracking regardless of outcome + if not paths: + return None + final_path = checkUniqueLogFile(out_path) + final_path.parent.mkdir(parents=True, exist_ok=True) + with open(final_path, 'w') as pipeline_log: + for path in paths: + section = path.parent.name # e.g. 'convert', 'index', 'search' + pipeline_log.write(f"\n{'=' * 80}\n# {section} ({path})\n{'=' * 80}\n\n") + if path.exists(): + pipeline_log.write(path.read_text()) + else: + pipeline_log.write(f'[log file missing: {path}]\n') + return final_path # -- _plainFormatter: internal helper to format log messages def _plainFormatter() -> PlainFormatter: diff --git a/src/comms/utils/modspec.py b/src/comms/utils/modspec.py new file mode 100644 index 0000000..65caa44 --- /dev/null +++ b/src/comms/utils/modspec.py @@ -0,0 +1,181 @@ +''' +comMS shared modification-spec and organism helpers + +Used by commands/config.py (persisting changes to a config file) and commands/index.py +(applying one-off, non-persisted overrides for a single `comms index` run). Every +function here is pure — it takes a value in and returns a new value out. +''' + +# -- Import external dependencies +import re + +# -- Import internal functions +from comms.utils.log import logMsg + +# -- Modification constants +CARBAMIDOMETHYL_MOD = 'C+57.0215' # static carbamidomethylation of Cys +MET_OX_MOD = '1M+15.9949' # variable Met oxidation +PHOSPHO_MOD = '1STY+79.966331' # variable STY phosphorylation +NCYC_MOD = '1Q-17.027' # N-terminal Gln cyclisation +NACE_MOD = '1X+42.011' # N-terminal protein acetylation +MANAGED_MOD_PATTERNS: dict[str, str] = { + r'^\d*C[+\-]': '--iodo / --no-iodo', + r'^\d*M\+15\.9949': '--ox / --no-ox', + r'^\d*STY\+79\.966331': '--phos / --no-phos', +} # mods that --custom is not allowed to duplicate (maps the residue/pattern that identifies each managed mod to its flag name) + +# -- Resolution constants +MZ_BIN_WIDTH_HIGH_RES = 0.02 # high-resolution instruments (default) +MZ_BIN_WIDTH_LOW_RES = 1.0005079 # low-resolution instruments +SCORE_FUNC_HIGH_RES = 'xcorr' # high-resolution instruments (default) +SCORE_FUNC_LOW_RES = 'combined-p-value' # low-resolution instruments + +# -- apply_mod: returns mods_spec string +def apply_mod(mods_spec: str, mod: str, exclusive_pattern: str | None = None) -> str: + ''' + Add or remove a mod entry in a Tide mods_spec string. + ''' + entries = [e.strip() for e in mods_spec.split(',') if e.strip()] + if exclusive_pattern: + pattern = re.compile(exclusive_pattern, re.IGNORECASE) + entries = [e for e in entries if not pattern.match(e)] + elif mod == '': + pass + else: + entries = [e for e in entries if e != mod] + if mod: + entries = [mod] + entries + return ','.join(entries) + +# -- apply_iodo: returns fixed_mods string +def apply_iodo(fixed_mods: str, iodo: bool) -> str: + ''' + Add or remove the carbamidomethylation Cys mod in a Tide fixed_mods string + ''' + entries = [e.strip() for e in fixed_mods.split(',') if e.strip()] + entries = [e for e in entries if e != CARBAMIDOMETHYL_MOD and e != 'C+0'] + if iodo: + entries = [CARBAMIDOMETHYL_MOD] + entries + else: + entries = ['C+0'] + entries # Crux automatically adds cysteine carbamidomethylation unless this string present + return ','.join(entries) + +# -- apply_custom_mod: returns custom_mods string +def apply_custom_mod(custom_mods: str, new_entry: str) -> str: + ''' + Add a custom mod entry to the custom_mods string, or clear all custom mods if new_entry is an empty string + ''' + if new_entry == '': + return '' + for pattern, flag_name in MANAGED_MOD_PATTERNS.items(): + if re.match(pattern, new_entry, re.IGNORECASE): + logMsg.warn(f'{new_entry} is managed by the {flag_name} flag, ignoring') + return custom_mods + entries = [e.strip() for e in custom_mods.split(',') if e.strip()] + if new_entry not in entries: + entries.append(new_entry) + out_str = ','.join(entries) + logMsg.debug(f'custom_mods updated: {out_str}') + return out_str + +# -- apply_organism: returns config dict with organism section replaced +def apply_organism(cfg: dict, organism: dict[str, str]) -> dict: + ''' + Replace the [organism] section of a config dict with the supplied dictionary. + ''' + cfg['organism'] = organism + logMsg.debug(f'organism section set to {organism}') + return cfg + +# -- parse_organism_arg: returns dict parsed from list of 'Key=Pattern' strings +def parse_organism_arg(pairs: list[str]) -> dict[str, str]: + ''' + Parse a list of 'Label=Pattern' strings into a dict. + ''' + result = {} + for item in pairs: + if '=' not in item: + logMsg.error(f'Invalid organism argument {item} (expected format: Organism=Pattern)') + raise SystemExit(1) + key, _, pattern = item.partition('=') + key = ''.join(key.split()) + pattern = ''.join(pattern.split()) + if not key: + logMsg.error(f'Empty label in organism argument: {item}') + raise SystemExit(1) + if not pattern: + logMsg.error(f'Empty pattern in organism argument: {item}') + raise SystemExit(1) + result[key] = pattern + return result + +# -- apply_protocol_flags: returns a config dict with protocol-level flags applied +def apply_protocol_flags( + cfg: dict, + *, + iodo: bool | None = None, + ox: bool | None = None, + phos: bool | None = None, + n_cyc: bool | None = None, + n_ace: bool | None = None, + clip_met: bool | None = None, + low_res: bool | None = None, + missed_cleavages: int | None = None, +) -> dict: + ''' + Apply protocol flags to a config dictionary and return it + iodo — owns the Cys slot in index.fixed_mods exclusively + ox — adds/removes 1M+15.9949 in index.mods_spec + phos — adds/removes 1STY+79.966331 in index.mods_spec + n_cyc — adds/removes 1Q-17.027 in index.nterm_peptide_mods_spec + n_ace — adds/removes 1X+42.011 in index.nterm_protein_mods_spec + low_res — sets search.mz_bin_width and search.score_function + missed_cleavages — sets index.missed_cleavages directly + ''' + cfg.setdefault('search', {}) + cfg.setdefault('index', {}) + cfg['index'].setdefault('fixed_mods', '') + cfg['index'].setdefault('nterm_peptide_mods_spec', '') + cfg['index'].setdefault('nterm_protein_mods_spec', '') + if iodo is not None: + cfg['index']['fixed_mods'] = apply_iodo(cfg['index'].get('fixed_mods', ''), iodo=iodo) + logMsg.debug(f'{"--iodo" if iodo else "--no-iodo"} applied: fixed_mods updated to {cfg["index"]["fixed_mods"]}') + if ox is not None: + spec = cfg['index'].get('mods_spec', '') + if ox: + cfg['index']['mods_spec'] = apply_mod(spec, mod=MET_OX_MOD) + else: + cfg['index']['mods_spec'] = apply_mod(spec, mod='', exclusive_pattern=r'^\d*M\+15\.9949') + logMsg.debug(f'{"--ox" if ox else "--no-ox"} applied: mods_spec updated to {cfg["index"]["mods_spec"]}') + if phos is not None: + spec = cfg['index'].get('mods_spec', '') + if phos: + cfg['index']['mods_spec'] = apply_mod(spec, mod=PHOSPHO_MOD) + else: + cfg['index']['mods_spec'] = apply_mod(spec, mod='', exclusive_pattern=r'^\d*STY\+79\.966331') + logMsg.debug(f'{"--phos" if phos else "--no-phos"} applied: mods_spec updated to {cfg["index"]["mods_spec"]}') + if n_cyc is not None: + spec = cfg['index'].get('nterm_peptide_mods_spec', '') + if n_cyc: + cfg['index']['nterm_peptide_mods_spec'] = apply_mod(spec, mod=NCYC_MOD) + else: + cfg['index']['nterm_peptide_mods_spec'] = apply_mod(spec, mod='', exclusive_pattern=r'^\d*Q\-17\.027') + logMsg.debug(f'{"--n-cyc" if n_cyc else "--no-n-cyc"} applied: nterm_peptide_mods_spec updated to {cfg["index"]["nterm_peptide_mods_spec"]}') + if n_ace is not None: + spec = cfg['index'].get('nterm_protein_mods_spec', '') + if n_ace: + cfg['index']['nterm_protein_mods_spec'] = apply_mod(spec, mod=NACE_MOD) + else: + cfg['index']['nterm_protein_mods_spec'] = apply_mod(spec, mod='', exclusive_pattern=r'^\d*X\+42\.011') + logMsg.debug(f'{"--n-ace" if n_ace else "--no-n-ace"} applied: nterm_protein_mods_spec updated to {cfg["index"]["nterm_protein_mods_spec"]}') + if low_res is not None: + cfg['search']['mz_bin_width'] = MZ_BIN_WIDTH_LOW_RES if low_res else MZ_BIN_WIDTH_HIGH_RES + cfg['search']['score_function'] = SCORE_FUNC_LOW_RES if low_res else SCORE_FUNC_HIGH_RES + logMsg.debug(f'{"--low-res" if low_res else "--high-res"} applied: mz_bin_width: {cfg["search"]["mz_bin_width"]}, score_function: {cfg["search"]["score_function"]}') + if clip_met is not None: + cfg['index']['clip_n_met'] = bool(clip_met) + logMsg.debug(f'{"--clip-met" if clip_met else "--no-clip-met"} applied: {cfg["index"]["clip_n_met"]}') + if missed_cleavages is not None: + cfg['index']['missed_cleavages'] = missed_cleavages + logMsg.debug(f'--missed-cleavages applied: {missed_cleavages}') + return cfg \ No newline at end of file diff --git a/src/comms/utils/readiness.py b/src/comms/utils/readiness.py index ba46841..9646612 100644 --- a/src/comms/utils/readiness.py +++ b/src/comms/utils/readiness.py @@ -12,17 +12,20 @@ def missing_requirements( has_sample_sheet, has_organism_prefix, multispecies, - has_organism_tags + has_organism_tags, + has_trfp, + has_crux, + has_r_deps, ) -> dict[str, list[str]]: '''Return, per command, the human-readable inputs still missing to run it.''' base = { - 'convert': [('data files', has_data)], - 'index': [('database', has_database)], - 'search': [('data files', has_data), ('database', has_database)], - 'rescore': [('database', has_database)] + [('organism patterns', has_organism_tags)], - 'lfq': [('sample sheet', has_sample_sheet), ('data files', has_data)], - 'quantify': [('database', has_database)], - 'report': [('sample sheet', has_sample_sheet), ('organism prefix', has_organism_prefix)], + 'convert': [('data files', has_data), ('ThermoRawFileParser', has_trfp)], + 'index': [('database', has_database), ('Crux', has_crux)], + 'search': [('data files', has_data), ('database', has_database), ('Crux', has_crux)], + 'rescore': ([('database', has_database), ('Crux', has_crux)] + ([('organism patterns', has_organism_tags)] if multispecies else [])), + 'lfq': [('sample sheet', has_sample_sheet), ('data files', has_data), ('Crux', has_crux)], + 'quantify': [('database', has_database), ('Crux', has_crux)], + 'report': [('sample sheet', has_sample_sheet), ('organism prefix', has_organism_prefix), ('R dependencies', has_r_deps)], } - base['pipeline'] = [item for items in base.values() for item in items] + base['pipeline'] = [item for cmd, items in base.items() if cmd != 'convert' for item in items] return {cmd: [name for name, ok in items if not ok] for cmd, items in base.items()} \ No newline at end of file diff --git a/src/comms/utils/settings.py b/src/comms/utils/settings.py index 113285f..ad45545 100644 --- a/src/comms/utils/settings.py +++ b/src/comms/utils/settings.py @@ -3,16 +3,19 @@ ''' # -- Import external dependencies -import tomllib +import tomllib, tomli_w from importlib.resources import files as pkg_files from pathlib import Path from platformdirs import user_config_dir from rich import print -from typing import Optional +from typing import Optional, TypeVar # Import internal classes/functions from comms.utils.log import logMsg +# -- Define TypeVar T +T = TypeVar('T') + # -- globalConfigPath: returns Path to OS-appropriate config file def globalConfigPath() -> Path: ''' @@ -51,6 +54,12 @@ def _loadTomlFile(path: Path) -> dict: with path.open('rb') as f: return tomllib.load(f) +# -- _writeConfigTo: writes a config dict to a given path +def _writeConfigTo(config: dict, path: Path) -> None: + path.parent.mkdir(parents=True, exist_ok=True) + with path.open('wb') as f: + tomli_w.dump(config, f) + # -- resolveConfig: returns (config, source) def resolveConfig(comms_dir: Optional[Path] = None) -> tuple[dict, str]: ''' @@ -68,6 +77,18 @@ def resolveConfig(comms_dir: Optional[Path] = None) -> tuple[dict, str]: logMsg.debug('Using bundled default config') return loadDefaultConfig(), 'bundled defaults' +# -- resolve_config_value: returns override if given, else the config.toml value at [section].key +def resolve_config_value(cfg: dict, section: str, key: str, override: Optional[T]) -> T: + ''' + Return override if given (not None), else cfg[section][key + ''' + if override is not None: + return override + try: + return cfg[section][key] + except KeyError: + raise KeyError(f'No value for [{section}].{key} in config, and no override given') from None + # -- initComms: returns None, but prints start-up message to terminal def initComms() -> None: ''' diff --git a/src/comms/utils/sheet.py b/src/comms/utils/sheet.py index 2cf4cad..fa49170 100644 --- a/src/comms/utils/sheet.py +++ b/src/comms/utils/sheet.py @@ -32,4 +32,26 @@ def render_sample_sheet(rows) -> str: lines.append('\t'.join( [r.sample_id, r.raw_file, r.treatment, r.fraction, replicate, r.batch] )) - return '\n'.join(lines) + '\n' \ No newline at end of file + return '\n'.join(lines) + '\n' + +# -- parse_sample_sheet: parse a TSV text into a list of SampleRow +def parse_sample_sheet(text: str) -> list['SampleRow']: + lines = [line for line in text.splitlines() if line.strip()] + if not lines: + return [] + header = [h.strip() for h in lines[0].split('\t')] + rows: list[SampleRow] = [] + for line in lines[1:]: + values = line.split('\t') + record = dict(zip(header, values)) + replicate_text = record.get('replicate', '').strip() + rows.append(SampleRow( + sample_id=record.get('sample_id', '').strip(), + raw_file=record.get('raw_file', '').strip(), + treatment=record.get('treatment', '').strip(), + fraction=record.get('fraction', '').strip(), + replicate=int(replicate_text) if replicate_text else None, + batch=record.get('batch', '').strip(), + replicate_overridden=bool(replicate_text), + )) + return rows \ No newline at end of file diff --git a/src/comms/utils/validate.py b/src/comms/utils/validate.py index a820b41..b4c050c 100644 --- a/src/comms/utils/validate.py +++ b/src/comms/utils/validate.py @@ -183,4 +183,14 @@ def _run(cmd: list[str]) -> Optional[str]: output = _run([mono, str(trfp_path), '--version']) if output is None: return None - return _parse_version(output) \ No newline at end of file + return _parse_version(output) + +# -- probe_crux: returns Path to Crux binary if found under bin_dir, else None (non-raising version for experiment GUI) +def probe_crux(bin_dir: Path) -> Optional[Path]: + result = _select_best(_find_all_crux(bin_dir), _get_crux_version) + return result[0] if result else None + +# -- probe_trfp: returns Path to ThermoRawFileParser binary if found under bin_dir, else None (non-raising version for experiment GUI) +def probe_trfp(bin_dir: Path) -> Optional[Path]: + result = _select_best(_find_all_trfp(bin_dir), _get_trfp_version) + return result[0] if result else None \ No newline at end of file diff --git a/tests/TESTS.md b/tests/TESTS.md index 1c901ed..427deb7 100644 --- a/tests/TESTS.md +++ b/tests/TESTS.md @@ -5,14 +5,16 @@ This document outlines the comMS test suite: its structure, shared fixtures, and - [Running the test suite](#running-the-test-suite) - [Test markers](#test-markers) - [Shared fixtures](#shared-fixtures) - - [Binary availability fixtures](#binary-availability-fixtures) + - [Bin directory and binary availability fixtures](#bin-directory-and-binary-availability-fixtures) - [Synthetic file fixtures](#synthetic-file-fixtures) - [Sample sheet fixtures](#sample-sheet-fixtures) - [Config fixtures](#config-fixtures) - [Synthetic Percolator results](#synthetic-percolator-results) - - [LFQ fixtures](#lfq-fixtures) - - [Rescore integration fixtures](#rescore-integration-fixtures) + - [PSM directory fixtures](#psm-directory-fixtures) - [GUI fixtures](#gui-fixtures) + - [Experiment context fixtures](#experiment-context-fixtures) + - [Experiment builder fixture](#experiment-builder-fixture) + - [Rescore/crux integration fixtures](#rescorecrux-integration-fixtures) - [Synthetic fixture generator](#synthetic-fixture-generator) - [Running standalone](#running-standalone) - [Synthetic proteome](#synthetic-proteome) @@ -21,21 +23,27 @@ This document outlines the comMS test suite: its structure, shared fixtures, and - [Real `.RAW` fixture](#real-raw--fixture) - [Unit tests](#unit-tests) - [`tests/unit/test_config.py`](#testsunittest_configpy) + - [`tests/unit/test_modspec.py`](#testsunittest_modspecpy) - [`tests/unit/test_context.py`](#testsunittest_contextpy) - [`tests/unit/test_experiment.py`](#testsunittest_experimentpy) - [`tests/unit/test_fasta.py`](#testsunittest_fastapy) + - [`tests/unit/test_installrdeps.py`](#testsunittest_installrdepspy) - [`tests/unit/test_lfq.py`](#testsunittest_lfqpy) + - [`tests/unit/test_license.py`](#testsunittest_licensepy) - [`tests/unit/test_parammedic.py`](#testsunittest_parammedicpy) - [`tests/unit/test_paths.py`](#testsunittest_pathspy) + - [`tests/unit/test_readiness.py`](#testsunittest_readinesspy) - [`tests/unit/test_report.py`](#testsunittest_reportpy) - [`tests/unit/test_rescore.py`](#testsunittest_rescorepy) - [`tests/unit/test_samples.py`](#testsunittest_samplespy) - [`tests/unit/test_settings.py`](#testsunittest_settingspy) + - [`tests/unit/test_sheet.py`](#testsunittest_sheetpy) - [`tests/unit/test_uninstall.py`](#testsunittest_uninstallpy) - [`tests/unit/test_validate.py`](#testsunittest_validatepy) - [`tests/unit/test_version.py`](#testsunittest_versionpy) - [`tests/unit/gui/test_gui_models.py`](#testsunitguitest_gui_modelspy) - [`tests/unit/gui/test_gui_panels.py`](#testsunitguitest_gui_panelspy) + - [`tests/unit/gui/test_gui_readiness.py`](#testsunitguitest_gui_readinesspy) - [`tests/unit/gui/test_gui_status.py`](#testsunitguitest_gui_statuspy) - [`tests/unit/gui/test_gui_widgets.py`](#testsunitguitest_gui_widgetspy) - [Integration tests](#integration-tests) @@ -44,6 +52,7 @@ This document outlines the comMS test suite: its structure, shared fixtures, and - [`tests/integration/test_pipeline.py`](#testsintegrationtest_pipelinepy) - [`tests/integration/test_trfp.py`](#testsintegrationtest_trfppy) - [R tests](#r-tests) + - [`tests/r/helper.R`](#testsrhelperr) - [`tests/r/test_utils_import.R`](#testsrtest_utils_importr) - [`tests/r/test_utils_normalise.R`](#testsrtest_utils_normaliser) @@ -83,8 +92,8 @@ Two custom markers control which tests are run depending on external binary avai | Marker | Requirement | |---|---| -| `crux` | Requires the Crux binary under `bin/` | -| `trfp` | Requires ThermoRawFileParser under `bin/` | +| `crux` | Requires the Crux binary under `tests/bin/` | +| `trfp` | Requires ThermoRawFileParser under `tests/bin/` | Tests decorated with these markers are skipped automatically (with a message) if the corresponding binary is not found, and the rest of the suite continue unaffected. To run only unmarked tests (i.e. those with no external dependency): @@ -98,16 +107,18 @@ uv run pytest -m "not crux and not trfp" ## Shared fixtures Shared fixtures are defined in `tests/conftest.py` and are described below. -### Binary availability fixtures -Two session-scoped fixtures are defined to locate external binaries, which both use the same regular expression/globbing logic as the source code. If a binary is absent, any test depending on that fixture will be skipped. +### Bin directory and binary availability fixtures +A session-scoped, autouse fixture (`_comms_bin_dir_env`) sets the `COMMS_BIN_DIR` environment variable to `tests/bin` for the whole test session, so every test — not just those requesting `crux_bin`/`trfp_exe` directly — resolves binaries against the bundled test `bin/` directory rather than a real installation. + +Two session-scoped fixtures locate external binaries within that directory, mirroring the resolution logic in `comms.utils.crux.findCrux` and `comms.utils.trfp.findTRFP`. If a binary is absent, any test depending on that fixture is skipped. Fixture | Description -- | -- -`crux_bin` | Resolves to Crux binary path under `bin/`; requires `pytest.mark.crux` -`trfp_exe` | Resolves to `ThermoRawFileParser.exe` under `bin`; requires `pytest.mark.trfp` +`crux_bin` | Resolves to the highest-versioned `crux*/bin/crux` under `tests/bin/`; requires `pytest.mark.crux` +`trfp_exe` | Resolves to `*/ThermoRawFileParser(.exe)` under `tests/bin/`; requires `pytest.mark.trfp` ### Synthetic file fixtures -Fixture | Description +Fixture | Description ---|--- `synthetic_fasta(tmp_path)` | Writes `synthetic_proteome.fasta` to a temporary directory and returns its path `synthetic_mzml(tmp_path)` | Writes `synthetic.mzML` to a temporary directory and returns its path @@ -116,55 +127,55 @@ Fixture | Description ### Sample sheet fixtures Fixture | Description ---|--- -`valid_sample_sheet(tmp_path)` | Two samples across two treatments with a `fraction` column (`WCL`) and optional `batch` column, one replicate each; written as TSV +`sample_sheet_factory(tmp_path)` | Returns a function `_make(fractions=('WCL',), batch=True)` that writes a minimal TSV sample sheet (two treatments, `MOCK`/`TREAT`, one replicate each, per given fraction) and returns its path. Replaces the old fixed `valid_sample_sheet*` fixtures, letting each test ask for exactly the fractions it needs `sample_sheet_missing_col(tmp_path)` | Missing the required `treatment` column; used to test validation errors `sample_sheet_duplicate_ids(tmp_path)` | Duplicate `sample_id` values; used to test duplicate detection ### Config fixtures Fixture | Description ---|--- -`isolated_config_dir(tmp_path, monkeypatch)` | Monkeypatches `globalConfigPath()` in both `settings` and `config` modules to point at a temporary directory, so tests don't access the real OS config file +`isolated_config_dir(tmp_path, monkeypatch)` | Monkeypatches `globalConfigPath()` in the `settings`, `config`, and `uninstall` modules to point at a temporary directory, so tests don't touch the real OS config file ### Synthetic Percolator results - Fixture | Description ---|--- `synthetic_percolator_results(tmp_path)` | Writes a minimal synthetic Percolator PSM file at `rescore/EUK/synthetic.EUK.percolator.target.psms.txt`, matching the per-organism subdirectory structure produced by `run_rescore` round 2 and bypassing the need to run Percolator on synthetic data (which does not provide enough PSMs for convergence) -### LFQ fixtures - -Fixture | Description ----|--- -`valid_sample_sheet_single_fraction(tmp_path)` | Single fraction (`WCL`), two treatments; used to test the single-fraction edge case in `_groupPsmsByFraction` -`valid_sample_sheet_multiple_fractions(tmp_path)` | Three fractions (`WCL`, `ECF`, `PUR`), two treatments, one replicate each; written as TSV -`single_fraction_psm_dir(tmp_path)` | Writes two synthetic Percolator PSM files (one fraction) to `comms/results/rescore/`, matching `valid_sample_sheet_single_fraction`; returns the directory path -`multi_fraction_psm_dir(tmp_path)` | Writes six synthetic Percolator PSM files (two per fraction) to `comms/results/rescore/`, matching `valid_sample_sheet_multiple_fractions`; returns the directory path - -### Rescore integration fixtures - -The following fixtures are defined within `tests/integration/test_pipeline.py` for use by the rescore integration test classes. - +### PSM directory fixtures Fixture | Description ---|--- -`two_organism_fasta(tmp_path)` | Writes a combined FASTA containing one TESTEUK protein, one TESTPRO protein, and one cRAP contaminant; returns the path -`synthetic_tide_search_dir(tmp_path)` | Writes a minimal synthetic Tide-search target PSM file to `tmp_path / 'search'` and returns the directory path, bypassing the need to run `tide-search` in rescore tests - -A module-level helper function `_write_per_organism_psm_files(rescore_dir, file_base, labels)` is also defined in `test_pipeline.py`. This is not a fixture but is used as a mock side effect within `TestRunRescoreMergedOutput` to write synthetic per-organism Percolator output files so that `_mergeRescoredPsms` has something to read without Percolator running. +`psm_dir_factory(tmp_path)` | Returns a function `_make(stems)` that writes one synthetic Percolator PSM file per given stem to `comms/results/rescore/` and returns that directory. Replaces the old fixed `single_fraction_psm_dir`/`multi_fraction_psm_dir` fixtures ### GUI fixtures -`QT_QPA_PLATFORM=offscreen` is set at the top of `conftest.py` so Qt widgets can be constructed and painted without a display +`QT_QPA_PLATFORM=minimal` is set at the top of `conftest.py` so Qt widgets can be constructed and painted without a display. A custom Qt message handler (installed at import time via `qInstallMessageHandler`) silently drops the small set of known-harmless warnings the offscreen/minimal platform plugin emits (e.g. "does not support grabbing the keyboard/mouse") and forwards everything else to stderr as Qt would by default. Fixture | Description -- | -- -`qapp` | Session-scoped `QApplication` (created with `QApplication.instance() or QApplication([])`) so GUI widgets can be built in tests; requested via `pytestmark = pytest.mark.usefixtures('qapp')` at the top of each GUI test module +`qapp` | Session-scoped `QApplication` (created with `QApplication.instance() or QApplication([''])`) so GUI widgets can be built in tests; requested via `pytestmark = pytest.mark.usefixtures('qapp')` at the top of each GUI test module ### Experiment context fixtures +Fixture | Description +-- | -- +`experiment_ctx(tmp_path, isolated_config_dir)` | Returns `ExperimentContext.resolve(tmp_path)`: an experiment context rooted at `tmp_path` with no `experiment.toml`, so config resolves to the bundled defaults and `bin_dir` is `None`. Also depends on `isolated_config_dir` so global-config fallback in tests never touches the real OS config -A bare experiment context is provided for the command-level integration tests, which now receive an ExperimentContext rather than a raw output directory. - +### Experiment builder fixture Fixture | Description -- | -- -`experiment_ctx(tmp_path)` | Returns ExperimentContext.resolve(tmp_path): an experiment context rooted at tmp_path with no experiment.toml, so config resolves to the bundled defaults and bin_dir is None +`experiment_builder(tmp_path, isolated_config_dir, sample_sheet_factory, psm_dir_factory)` | Returns a chainable builder for composing a `comms/` directory from only the pieces a test needs: `.with_sample_sheet(fractions, batch)` writes a sample sheet and records it under `[files]`; `.with_stage_output(stage, files)` writes placeholder output for a pipeline stage (PSM files via `psm_dir_factory` for `'rescore'`, otherwise empty placeholder files); `.with_metadata(**kwargs)` merges arbitrary sections into `experiment.toml`; `.build()` writes `experiment.toml` and returns `(root, ctx)` via `ExperimentContext.resolve` + +### Rescore/crux integration fixtures +The following fixtures are defined within `tests/integration/test_pipeline.py` and `tests/integration/test_crux.py` for use by their respective integration test classes. + +Fixture | Description +---|--- +`two_organism_fasta(tmp_path)` | Writes a combined FASTA containing one TESTEUK protein, one TESTPRO protein, and one cRAP contaminant; returns the path +`synthetic_tide_search_dir(tmp_path)` | Writes a minimal synthetic Tide-search target PSM file to `tmp_path / 'search'` and returns the directory path, bypassing the need to run `tide-search` in rescore tests +`pipeline_index(crux_bin, tmp_path_factory)` | Module-scoped: builds one shared Tide index for all `TestRunSearch*` tests in `test_pipeline.py`, skipping the module's tests if `run_index` fails +`pipeline_search(crux_bin, pipeline_index, tmp_path_factory)` | Module-scoped: runs `run_search` once and shares `(search_dir, fasta, work)` across `TestRunRescore` tests, skipping if `run_search` fails +`built_index(crux_bin, tmp_path_factory)` | Module-scoped, in `test_crux.py`: builds one shared Tide index via `tideIndex` for all `TestTideSearch` tests +`search_results(crux_bin, built_index, tmp_path_factory)` | Module-scoped, in `test_crux.py`: runs `tideSearch` once and returns `(out_dir, target_file, fasta)` + +Module-level helper functions `_write_combined_percolator_output`, `_write_split_psm_files`, and `_write_combined_psm` are also defined in `test_pipeline.py`/`test_rescore.py`. These are not fixtures but are used as mock side effects to write synthetic per-round output so that downstream logic (`_splitPsmsByOrganism`, the per-organism Percolator round) has something to read without Percolator or Tide-search actually running. ---

^ Back to top

@@ -182,7 +193,7 @@ By default, the script will write files to the directory containing the script. python tests/fixtures/generate_fixtures.py path/to/output/dir ``` ### Synthetic proteome -The synthetic proteome is written to `synthetic_proteome.fasta`. It contains five protein sequences, each with one or two tryptic peptides which exclusively map to the source protein. The protein IDs and peptide sequences are: +The synthetic proteome is written to `synthetic_proteome.fasta`. It contains five protein sequences; PROT1's sequence (`ACDEFGHIKLMNPQRSTVWYK`) is a single tryptic run yielding two of the target peptides, so that a multi-peptide protein is represented without duplicating sequence content. The protein IDs and peptide sequences are: Protein ID | Tryptic peptides | ---|--- @@ -222,7 +233,7 @@ Peptide mass = sum(residues) + water b-ions and y-ions are singly charged and skip the terminal ions (b1 and y1), following the standard convention. ### Real `.RAW ` fixture -Valid synthetic `.RAW` file cannot be generated without the ThermoFisher vendor SDK, therefore integration tests which would require a `.RAW` file are gated behind the `REAL_RAW_FIXTURE` guard. To run these tests, place a valid file at: +Valid synthetic `.RAW` file cannot be generated without the ThermoFisher vendor SDK, therefore integration tests which would require a `.RAW` file are gated behind the `REAL_RAW_FIXTURE` guard (a constant defined in `test_trfp.py` and imported by `test_convert.py`). To run these tests, place a valid file at: ``` tests/fixtures/real_sample.RAW ``` @@ -235,47 +246,53 @@ If this file is absent, the relevant tests are skipped automatically. Unit tests cover logic in isolation, i.e. they do not require external binaries and do not write to the local filesystem beyond `tmp_path`. All config-touching tests use the `isolated_config_dir` fixture described [above](#config-fixtures). ### `tests/unit/test_config.py` -Unit tests covering `src/comms/commands/config.py`, and indirectly `src/comms/utils/settings.py`: +Unit tests covering `src/comms/commands/config.py` following the config system overhaul (mod/protocol-flag logic has moved to `comms/utils/modspec.py`, covered separately below): Class | Test description -- | -- -`TestLoadDefaultConfig` | returns a dict; contains expected top-level sections; search section has required keys; `fragment_tolerance_da` key has been removed; key values have correct types; default score function is xcorr; default mz_bin_width is high-res; default fixed_mods does not contain carbamidomethyl -`TestFlatten` | flat dict unchanged; nested dict flattened; deeply nested; mixed depth; empty dict; default config flattens without error -`TestConfigCheck` | returns `True`/`False` correctly for exists/absent file under both `exists=True` and `exists=False` modes -`TestWriteLoadConfig` | round-trip preserves content; raises `FileNotFoundError` when no config present -`TestApplyMod` | adds mod to empty spec; adds to existing spec; prepends mod; duplicate entry not added; removal with exclusive pattern; exclusive pattern replaces on add; no leading/trailing commas; no double commas; empty mod with no pattern is no-op; removal of absent mod is no-op -`TestApplyIodo` | adds carbamidomethyl to empty and non-empty `fixed_mods` string; prepends; removes carbamidomethyl; no-op when not present; idempotent; result has no count prefix; no leading/trailing commas; no double commas -`TestApplyCustom` | adds entry to empty string; adds to existing; empty string clears all; duplicate not added; managed Met/Cys/phos mods rejected with warning; unmanaged entry accepted; no leading/trailing commas -`TestApplyOrganism` | sets organism section; replaces existing; does not touch other sections; empty dict clears; returns cfg -`TestApplyProtocolFlags` | iodo/low_res/mbr None leaves keys unchanged; low_res True/False sets bin width and score function; combined iodo+low_res; only relevant keys touched; non-search sections untouched -`TestApplyProtocolFlagsMods` | iodo True/False writes to `fixed_mods`, not `mods_spec`; ox/phos/n_cyc/n_ace True adds correct mod to correct key; False removes; None is no-op; n_cyc/n_ace do not touch mods_spec; ox and iodo coexist across different keys; iodo does not remove ox; all flags None changes nothing -`TestParseOrganismArg` | single and multiple pairs; strips whitespace; preserves regex chars; raises `SystemExit` on no `=`, empty key, empty pattern; returns dict; empty list returns empty dict; `=` in pattern preserved -`TestConfigInit` | creates config file; file is valid TOML; does not overwrite existing -`TestConfigExists` | exits nonzero when absent; does not raise when present -`TestConfigVerify` | valid config passes; missing key exits nonzero; exits nonzero when no config -`TestConfigReset` | `--force` restores defaults; without force prompts; confirms and resets on accept -`TestConfigSet` | creates config if absent; created config is valid TOML; all named mod flags add/remove correct mod in correct key; idempotent for ox/phos/n_cyc/n_ace; n_cyc/n_ace do not change mods_spec; custom adds entry; custom is additive; custom empty string clears; custom managed mod not added; iodo flags unchanged from original tests; low_res/organism/mbr unchanged from original tests; combined flags work together; no-flags exits nonzero; all unrelated config keys unchanged after set -`TestResolveConfigTarget` | None resolves to the global user config path; "global"/"GLOBAL" resolve case-insensitively to the global path; any other string is returned verbatim as a Path -`TestConfigSetLocalTarget` | writes to the supplied local path and produces valid TOML; the global user config is left untouched when a local target is given +`TestFlatten` | flat dict unchanged; nested dict flattened; deeply nested; empty dict; the bundled default config flattens without error +`TestPrintTable` | renders without raising when current matches defaults; renders without raising when a value diverges from defaults +`TestPrintDiffSummary` | no-changes case prints a "No changes" message; a changed key shows both old and new values; multiple changes are all reported (one ✓ per change) +`TestResolveOrCreate` | `use_global=True` resolves to the global config path; `use_global=False` resolves to `/comms/config.toml`; creates a defaults-derived file if none exists; raises `SystemExit` when both a bare `config.toml` and a nested `comms/config.toml` are present (ambiguous); prefers the bare config when only that is present; prompts and creates the nested config when neither is present (accepted via `_confirm`); declining the prompt exits with code 0 +`TestConfigList` | prints the config path and a table (`config_list`) +`TestConfigVerify` | a config created from defaults passes; missing a required key exits non-zero; an unexpected/extra key exits non-zero +`TestConfigReset` | `--force` resets to bundled defaults without prompting; without force, prompts via `_confirm` and exits with code 0 only when declined +`TestConfigSet` | all-`None` flags returns `False` and writes nothing; a protocol flag (`iodo`) round-trips into `index.fixed_mods`; `organism` round-trips into the `organism` section; `custom` round-trips into `index.custom_mods`; direct flags (`gzip`, `threads`, `picked_protein`, `measure`, `lfc_threshold`) round-trip into their respective sections (parametrised); unrelated sections (`search`, `rescore`, `quantify`, `convert`) are untouched by an unrelated flag; a diff summary (✓) is printed after a successful set +`TestConfigSetLocalTarget` | writes to the local `/comms/config.toml`, not the global user config, which remains untouched + +--- +### `tests/unit/test_modspec.py` +Unit tests covering `src/comms/utils/modspec.py` — the modification-spec and protocol-flag logic previously tested as part of `test_config.py`, now in its own module alongside the mod-name constants (`CARBAMIDOMETHYL_MOD`, `MET_OX_MOD`, `PHOSPHO_MOD`, `NCYC_MOD`, `NACE_MOD`) and resolution constants (`MZ_BIN_WIDTH_HIGH_RES`/`LOW_RES`, `SCORE_FUNC_HIGH_RES`/`LOW_RES`): + +Class | Test description +-- | -- +`TestApplyMod` | adds mod to empty spec; adds to existing spec; prepends mod; duplicate entry not added; removal via exclusive pattern; exclusive pattern replaces on add; no leading/trailing commas; removal of an absent mod is a no-op +`TestApplyIodo` | adds carbamidomethyl when `iodo=True`; removes it when `iodo=False`; `iodo=False` on an empty spec adds `C+0`; idempotent (re-applying does not duplicate) +`TestApplyCustomMod` | adds entry to empty string; empty string clears all; duplicate not added; managed mods (Met ox, carbamidomethyl, STY phospho) are rejected and silently dropped (parametrised); `MANAGED_MOD_PATTERNS` is a `dict[str, str]` of pattern to flag name +`TestApplyOrganism` | sets the `organism` section; replaces an existing one; does not touch other sections +`TestParseOrganismArg` | single and multiple `key=value` pairs; strips whitespace; raises `SystemExit` on a missing `=`, empty key, or empty pattern; empty list returns an empty dict +`TestApplyProtocolFlags` | `low_res=True`/`False` sets `mz_bin_width` and `score_function` to the low-/high-res constants; `None` flags are a no-op; **new:** `clip_met=True`/`False` sets `index.clip_n_met` as a genuine bool (never a string) and `None` is a no-op; **new:** `missed_cleavages` sets `index.missed_cleavages` and `None` is a no-op; `clip_met` and `missed_cleavages` coexist and don't disturb mod-flag handling (`ox` still writes to `mods_spec`) --- ### `tests/unit/test_context.py` -Unit tests covering `src/comms/utils/context.py`. +Unit tests covering `src/comms/utils/context.py`: Class | Test description -- | -- `TestNormaliseDirs` | a plain root returns `(root, root/comms)`; a path whose final component is `comms` returns `(parent, comms)`; a directory containing `experiment.toml` directly is treated as a `comms` directory and returns `(parent, dir)` `TestResolve` | with no `experiment.toml` present, config falls back to the bundled default, `bin_dir` is `None`, and root equals the input directory; a local `comms/config.toml` is preferred and `config_source` begins with `"local"`; a `bin_dir` set in `experiment.toml` is parsed to a `Path` and exposed on the context -`TestExperimentContextProperties` | each stored-input `@property` (`data_files`, `database`, `sample_sheet`, `analysis_mode`, `multispecies`, `organism_prefix`, `ref_info`, `cont_csv`) returns the correct typed value when the corresponding `metadata` key is present, and `None` / empty list when absent; properties are **not** constructor parameters — they are read from `self.metadata` -`TestResultsDir` | returns the canonical `/comms/results/` path; the path is correct even when the directory does not yet exist -`TestChoose` | override is returned when given; a warning containing the label is logged when override differs from stored; no warning when override equals stored; stored is returned when no override; raises `SystemExit` when neither is supplied; raises `SystemExit` when resolved path is missing and `must_exist=True`; returns a non-existent path without raising when `must_exist=False` +`TestResolveResultsInput` | delegates to `results_dir` for the default (no-override) case; `must_exist=False` suppresses the existence check even though the directory doesn't exist; an explicit override wins over the computed path +`TestResolveReport` | `report_enabled=True` and no override means do-not-skip (`resolve_report` returns `False`); `report_enabled=False` means skip (`True`); `report_enabled=None` (unset) defaults to do-not-skip; an override agreeing with the context logs no warning; an override disagreeing with the context wins and logs a warning; the polarity is pinned explicitly — `resolve_report` returns *skip*, not *enabled* +`TestExperimentContextProperties` | each stored-input `@property` (`data_files`, `database`, `sample_sheet`, `analysis_mode`, `multispecies`, `organism_prefix`, `ref_info`, `cont_csv`, `report_enabled`) returns the correct typed value when the corresponding metadata key is present, and `None`/empty list when absent; `multispecies` falls back to checking whether `config['organism']` is non-empty when no `[experiment].analysis` mode is stored; properties are read from `self.metadata`, not constructor parameters +`TestResultsDir` | `results_dir(ctx, command)` returns the canonical `/comms/results/` path; correct even when the directory does not yet exist +`TestChoose` | `_choose` returns the override when given; logs a warning naming the label when override differs from stored; no warning when override equals stored; returns stored when no override; raises `SystemExit` when neither is supplied; raises `SystemExit` when the resolved path is missing and `must_exist=True`; returns a non-existent path without raising when `must_exist=False` `TestCheckFiles` | all-existing files returns a `list[Path]` equal to the inputs (not `True`); a missing file raises `SystemExit`; empty input returns an empty list (not `True`) `TestResolveDataFiles` | stored list returned when no override; override wins and a warning is logged when it differs from stored; raises `SystemExit` when neither stored nor override is present; raises `SystemExit` when any file is missing; result is a `list[Path]` -`TestResolveMzmlFiles` | explicit override list returned directly; raises when override file missing; globs `*.mzML` and `*.mzML.gz` from the convert results directory when no override; raises when no files found in convert directory -`TestResolveSingleFileInputs` | `resolve_database` and `resolve_sample_sheet` return stored value; override wins with warning; raises when neither supplied; raises when resolved file is missing -`TestResolveOrganismPrefix` | stored prefix returned; override wins with warning; raises when neither is available; return type is `str` +`TestResolveMzmlFiles` | explicit override list returned directly; raises when an override file is missing; globs `*.mzML` and `*.mzML.gz` from the convert results directory when no override; raises when no files found in the convert directory +`TestResolveSingleFileInputs` | `resolve_database` and `resolve_sample_sheet` return the stored value; override wins with a warning; raises when neither is supplied; raises when the resolved file is missing +`TestResolveOrganismPrefix` | stored prefix returned; override wins with a warning; raises when neither is available; return type is `str` --- @@ -286,7 +303,7 @@ Class | Test description -- | -- `TestLaunchExperimentGui` | exits with `run_app`'s return code; logs a launch message; `logMsg` instance is named `'experiment'` `TestMainWindowCloseLogging` | closing the window logs a message containing "closed"; `logMsg` instance is named `'experiment'` -`TestRunExperimentHeadless` | with `typer.prompt`/`typer.confirm` patched and a temporary `.mzML` file, writes `sample_sheet.tsv`, `config.toml` and `experiment.toml` under `/comms/`; the prompt sequence now includes the database FASTA prompt (between bin-dir and treatments) and an explicit data-file list via `_prompt_list('data file')` (between the input directory and per-file assignment); requires at least one treatment and one fraction (exits non-zero otherwise); records a `bin_dir` in `experiment.toml` only when one is supplied +`TestRunExperimentHeadless` | with `typer.prompt`/`typer.confirm` patched and a temporary `.mzML` file, writes `sample_sheet.tsv`, `config.toml` and `experiment.toml` under `/comms/`; the prompt sequence includes the combined database FASTA prompt (between bin-dir and treatments) and an explicit data-file list via `_prompt_list('data file')` (between the input directory and per-file assignment); records a `bin_dir` in `experiment.toml` only when one is supplied --- @@ -302,6 +319,18 @@ Class | Test description --- +### `tests/unit/test_installrdeps.py` +Unit tests covering `src/comms/utils/installrdeps.py` — checking, printing and installing the R package dependencies used by the `report` command. No R installation is required; `shutil.which`, `subprocess.run`/`Popen` are mocked throughout. + +Class | Test description +-- | -- +`TestCheckRDependencies` | returns `None` when `Rscript` is not on `PATH`; parses a well-formed JSON `{installed, missing}` payload from stdout; a non-zero return code returns `None`; malformed JSON returns `None` +`TestPrintDependencyTable` | only an "Installed" line is logged when nothing is missing; both "Installed" and "Missing" lines (plus the `r-utils install` hint) are logged when something is missing; nothing is logged when both lists are empty +`TestInstallRDependencies` | returns `False` when `Rscript` is not on `PATH`; streams `Popen` stdout lines to the log as they arrive and returns `True` on a zero exit code; a non-zero exit code returns `False` +`TestInstallRDependenciesTerminal` | nothing missing informs the user without prompting; something missing and the user confirms (`logMsg.input` returns `'y'`) calls `install_r_dependencies`; the user declining (`'n'`) does not call it; malformed check output logs an error without raising + +--- + ### `tests/unit/test_lfq.py` Unit tests covering `_groupPsmsByFraction` in `src/comms/commands/lfq.py`. No external binaries are required. @@ -311,6 +340,15 @@ Class | Test description --- +### `tests/unit/test_license.py` +Unit tests covering `src/comms/commands/license.py`: + +Class | Test description +-- | -- +`TestPrintLicense` | raises `SystemExit` with code 0; reads and pages the license file without raising (`pydoc.pager` is mocked and asserted called once) + +--- + ### `tests/unit/test_parammedic.py` Unit tests covering `_parseParamMedicOutput` and `_runParamMedic` in `src/comms/commands/search.py`. No external binaries are required; `cruxutil.paramMedic` is mocked throughout `TestRunParamMedic`. @@ -326,61 +364,85 @@ Unit tests covering `src/comms/utils/paths.py`: Class | Test description -- | -- -`TestGenerateOutputFileStructure` | creates the expected `comms/results//` subdirectory; creates directories if absent; returns existing path unchanged if already correct; works for all supported commands -`TestCheckUniqueFileName` | returns expected base name when no conflict; increments suffix on conflict; increments correctly through multiple conflicts; correct naming patterns for all commands (`search`, `quantify`, `rescore`, `report`); returned path is within `out_dir` +`TestGenerateOutputFileStructure` | creates the expected `comms/results//` subdirectory; creates directories if absent; returns existing path unchanged if already correct; works for `convert`, `index`, `search`, `rescore`, `quantify` (parametrised) +`TestCheckUniqueFileName` | returns the expected base name when there's no conflict; increments the numeric suffix on conflict, and again through multiple conflicts; `quantify` naming (`.spectral-counts.txt`); `rescore` naming (`.percolator.psms.txt`); `report` naming takes a `fmt` kwarg instead of `orig_name` (`comms-report.`); an unrecognised command still produces a usable fallback name (`comms--output...`) rather than raising; the returned path's parent is `out_dir` `TestRepoBinDir` | an explicit `experiment_bin_dir` takes precedence over everything and is returned unchanged; explicit value beats `COMMS_BIN_DIR` even when both are set; `COMMS_BIN_DIR` is used when no explicit value is given; falls back to repo-root `bin/` path when neither is set; `experiment_bin_dir=None` is equivalent to omitting the argument; returns a `Path` object --- -### `tests/unit/test_report.py` -Unit tests covering `src/comms/commands/report.py`: +### `tests/unit/test_readiness.py` +Unit tests covering `src/comms/utils/readiness.py` — the command-readiness gap calculator that backs the GUI readiness panel (and, indirectly, the CLI's pre-flight checks): Class | Test description -- | -- -`TestResolveRScript` | returns a `Path`; path ends with the requested script name; auxiliary scripts under `aux/` subdirectory are resolved correctly -`TestWriteIndex` | creates `index.md`; contains all section names; failed sections marked FAILED; passed sections marked ✓; parameters block is included -`TestRunReportValidation` | raises `SystemExit` when no spectral-counts files in quantify directory; raises `SystemExit` when output directory exists without `--overwrite`; raises `SystemExit` when Rscript binary is not on PATH; silently drops concordance when `--lfq-dir` absent; creates output directory; writes `index.md`; passes `lfc_threshold` and `fdr_threshold` as positional args to the `da` section; `logMsg` instance is named `'report'` +`TestMissingRequirements` | `missing_requirements(**all_true)` returns an empty gap list for every command in `COMMANDS`; a missing `has_data` flag surfaces "data files" only in `convert`/`search`/`lfq`, not `index`/`quantify`; a missing `has_crux` flag surfaces "Crux" in every Crux-dependent command (`index`, `search`, `rescore`, `lfq`, `quantify`) but not `convert`/`report`; a missing `has_trfp` flag surfaces "ThermoRawFileParser" only in `convert`; a missing `has_r_deps` flag surfaces "R dependencies" only in `report`; `pipeline`'s gap list is exactly the union of every other command's gaps; rescore requires "organism patterns" only when `multispecies=True` (pinned explicitly as a regression guard for the §1.2 assumption) -N.B. `_run_r_section` and `shutil.which` are mocked throughout — no R installation is required. +--- + +### `tests/unit/test_report.py` +Unit tests covering helper functions and `run_report` in `src/comms/commands/report.py`. This module now tracks and reports per-organism outcomes for each report section, not just a single pass/fail per section. + +Class | Test description +-- | -- +`TestResolveRScript` | returns a `Path`; path ends with the requested script name; auxiliary scripts under `aux/` (e.g. `aux/ev-markers.R`) are resolved correctly +`TestWriteIndex` | creates `index.md`; contains all section names; a failed section is marked `FAILED`; a passed section is marked with `✓`; a `partial` section (mixed per-organism outcomes) is marked `PARTIAL`; a `skipped` section is marked `SKIPPED`; per-organism outcome lines are nested under their parent section line in the rendered output; the parameters block is included +`TestRunReportValidation` | raises `SystemExit` when no spectral-counts files are in the quantify directory; raises `SystemExit` when the output directory exists without `--overwrite`; raises `SystemExit` when the `Rscript` binary is not on `PATH`; the `concordance` section is silently dropped when no `--lfq-dir` is available (falls back to the conventional `comms/results/lfq` location, and only drops the section if that's also absent); `ref_info` falls back to the value stored on the experiment context when not passed explicitly; creates the output directory and writes `index.md`; passes `lfc_threshold` and `fdr_threshold` as positional args (not kwargs) to the `da` section; `logMsg` instance is named `'report'`; a config override (e.g. `lfc_threshold`) writes a `report.config.toml` sidecar recording the overridden value, while an unmodified run (all overrides `None`) writes no sidecar +`TestReadStatus` | a missing `_status.json` returns `({}, {})`; malformed JSON returns `({}, {})`; a well-formed file returns its `organisms` and `reasons` dicts; unrecognised status values are dropped from the parsed `organisms` dict +`TestSectionStatus` | no organisms and `proc_ok=True` is `'skipped'`; no organisms and `proc_ok=False` is `'failed'`; all-`ok` organisms is `'succeeded'`; a mix of `ok`/`failed` is `'partial'`; all-`failed` is `'failed'`; all-`skipped` (none ok or failed) is `'skipped'` +`TestLogOrganismOutcomes` | logs one line per organism; includes the failure reason in parentheses when present; omits the empty parenthetical when there's no reason --- ### `tests/unit/test_rescore.py` -Unit tests covering helper functions in `src/comms/commands/rescore.py`. No external binaries are required; tests use synthetic PSM files written directly to `tmp_path`. +Unit tests covering helper functions in `src/comms/commands/rescore.py`. No external binaries are required; tests use synthetic PSM files written directly to `tmp_path`. Note the organism split now happens on **Tide-search** target/decoy output (`_splitPsmsByOrganism`, ahead of the per-organism Percolator round), not on already-rescored Percolator PSMs — the protein-ID column position is located dynamically via `_findProteinIdsIndex` rather than assumed to be the last column. Class | Test description -- | -- -`TestParseOrganismTags` | parses two-organism comma-separated string; parses single-organism string; strips internal and leading/trailing whitespace; preserves regex characters in values; raises `SystemExit` on odd item count, single item, or empty string; returns `dict[str, str]`; keys and values are strings -`TestClassifyPsmRow` | returns a `list` for matching rows; returns `['EUK']` for a matching EUK row; returns `['PRO']` for a matching PRO row; returns the string `'contaminants'` for an unmatched row; returns `'contaminants'` for an empty row; uses the last tab-delimited column as the protein ID; returns a list with multiple labels when the protein ID matches more than one organism tag -`TestSplitPsmsByOrganism` | returns `True` on success; creates per-organism target files in labelled subdirectories; creates per-organism decoy files; EUK file contains only EUK rows; PRO file contains only PRO rows; contaminant rows go to a `contaminants/` bucket; header is preserved in each output file; returns a bool without raising when the target file is missing; skips a missing decoy file gracefully and still succeeds for the target; output files are non-empty; with `shared_policy='drop'`, rows matching more than one organism are excluded from all output files; with `shared_policy='include'`, rows matching more than one organism appear in all matching output files +`TestParseOrganismTags` | parses two-organism comma-separated string; parses single-organism string; strips internal and leading/trailing whitespace; preserves regex characters in values; raises `SystemExit` on odd item count, single item, or empty string; returns `dict[str, str]` +`TestClassifyPsmRow` | `_classifyPsmRow(row, id_index, organism_tags)` returns a `list` for matching rows; returns `['EUK']`/`['PRO']` for matching rows; returns `['contaminants']` for an unmatched or empty row; uses the column at the supplied `id_index` for the protein ID (rather than assuming the last column); returns a list with multiple labels when the protein ID matches more than one organism tag +`TestFindProteinIdsIndex` | `_findProteinIdsIndex(header, 'protein id')` returns the correct 0-based column index for a mid-row column; returns `0` when it's the first column +`TestSplitPsmsByOrganism` | operates on a combined Tide-search target/decoy file pair; returns `True` on success; creates per-organism target and decoy files under labelled subdirectories, named `.