diff --git a/.gitignore b/.gitignore index fb4caa0..19cce1a 100644 --- a/.gitignore +++ b/.gitignore @@ -10,6 +10,7 @@ src/xchemalign.egg-info/ src/xchemalign/_version.py test-data/outputs test-data/outputs/upload_* +!test-data/processing_stats/*.log *.ipynb venv* .vscode diff --git a/config/pdb-depo/xia2-multiplex.cif b/config/pdb-depo/xia2-multiplex.cif new file mode 100644 index 0000000..cc91e7e --- /dev/null +++ b/config/pdb-depo/xia2-multiplex.cif @@ -0,0 +1,22 @@ +data_software +# +loop_ +_software.pdbx_ordinal +_software.name +_software.classification +_software.type +_software.version +_software.date +_software.location +_software.description +_software.pdbx_reference_DOI +1 xia2 'data reduction' package ? ? https://github.com/xia2/xia2 'An expert system to perform automated X-ray diffraction data processing' 10.1107/S0021889809045701 +2 DIALS 'data reduction' package ? ? https://github.com/dials/dials 'Data processing and integration within the DIALS software package' 10.1107/S2059798317017235 +3 xia2 'data scaling' package ? ? https://github.com/xia2/xia2 'An expert system to perform automated X-ray diffraction data processing' 10.1107/S0021889809045701 +4 DIALS 'data scaling' package ? ? https://github.com/dials/dials 'Data processing and integration within the DIALS software package' 10.1107/S2059798317017235 +5 DIMPLE phasing program ? ? https://github.com/ccp4/dimple 'Automated macromolecular crystallography pipeline for rapid refinement and ligand screening.' ? +6 REFMAC phasing package ? ? https://www.ccp4.ac.uk/ 'An expert system to perform automated X-ray diffraction data processing' 10.1107/S0907444996012255 +7 DIMPLE 'model building' program ? ? https://github.com/ccp4/dimple 'Automated macromolecular crystallography pipeline for rapid refinement and ligand screening.' ? +8 Coot 'model building' package ? ? https://github.com/pemsley/coot 'Coot is for macromolecular model building, model completion and validation, particularly suitable for protein modelling using X-ray data.' 10.1107/S0907444996012255 +9 gemmi 'data extraction' library ? ? https://github.com/project-gemmi/gemmi 'A Python library developed primarily for use in the field of macromolecular crystallography (MX)' 10.21105/joss.04200 +# diff --git a/docs/source/DEV-GUIDE.md b/docs/source/DEV-GUIDE.md index 881d72a..85c7e9e 100644 --- a/docs/source/DEV-GUIDE.md +++ b/docs/source/DEV-GUIDE.md @@ -59,10 +59,10 @@ There is an environment at Diamond where users run the XChem align code on their This can be found on the Diamond file system at `/dls/science/groups/i04-1/software/xchem-align`. To roll out a new version of this: -1. Check that the repos are up to date on the `master` (XCA) and `main` (LNA) branches. +1. Check that the repo is up to date on the `master` branch. 2. Test locally 3. Tag the XCA repo and push the tag: `git tag 1.2.3` and `git push origin 1.2.3` (using the appropriate tag number) -4. ssh to Diamond and move into the `/dls/science/groups/i04-1/software/xchem-align` dir +4. ssh to Diamond. Then ssh to wilson. Then run `srun --nodes=1 --ntasks=4 --partition=cs05r --pty bash` and move into the `/dls/science/groups/i04-1/software/xchem-align` dir 5. `git pull` - update the repo 6. `git tag` - check you have the expected tag 7. `rm -rf env_xchem_align` - remove the old conda environment diff --git a/pyproject.toml b/pyproject.toml index 0094c60..6f1bbf2 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -28,6 +28,7 @@ dependencies = [ 'requests-toolbelt == 1.0.0', # for pdb deposition 'mmcif-gen == 1.2.2', + 'beautifulsoup4 == 4.14.3', # from lna 'loguru == 0.7.3', 'pydantic==2.6.0', diff --git a/src/pdbdepo/pdb_deposition.py b/src/pdbdepo/pdb_deposition.py index d36067f..22b8593 100644 --- a/src/pdbdepo/pdb_deposition.py +++ b/src/pdbdepo/pdb_deposition.py @@ -350,6 +350,10 @@ def process_input( data_processing_prog = row.get(Constants.SOAKDB_COL_DATA_PROCESSING_PROGRAM) if data_processing_prog == 'dials': data_processing_prog = 'xia2-dials' + elif data_processing_prog and data_processing_prog.startswith('xia2.multiplex'): + # value seen in the wild is "xia2.multiplex sample_group" (a free-text suffix + # identifying the multi-crystal grouping is appended after the program name) + data_processing_prog = 'xia2-multiplex' if data_processing_prog and data_processing_prog != 'None': data_processing_prog = data_processing_prog.lower() info('data processing was done with ' + data_processing_prog) diff --git a/src/pdbdepo/scrape_processing_stats.py b/src/pdbdepo/scrape_processing_stats.py index 11e0302..0f64293 100644 --- a/src/pdbdepo/scrape_processing_stats.py +++ b/src/pdbdepo/scrape_processing_stats.py @@ -2,6 +2,7 @@ import re from pathlib import Path import gemmi +from bs4 import BeautifulSoup from gemmi import cif from xchemalign import utils @@ -176,6 +177,98 @@ def handle_xia_3dii(file): return handle_text_file(file, d_xia_3dii) +# maps the row label in the xia2.multiplex HTML report's "Overall" table to the reflns key. +# the "Resolution" and "Completeness" rows need bespoke handling (range splitting / '%' stripping) +# so are not included here. +_XIA2_MULTIPLEX_ROW_KEYS = { + 'Observations': KEY_REFLNS_NUM_MEASURED, + 'Unique reflections': KEY_REFLNS_NUM_OBSERVED, + 'Multiplicity': KEY_REFLNS_PDBX_REDUNDANCY, + 'Mean I/σ(I)': KEY_REFLNS_NETI_OVER_SIGMA, + 'Rmerge': KEY_REFLNS_PDBX_RMERGE_I_OBS, + 'Rmeas': KEY_REFLNS_PDBX_RRIM_I_ALL, + 'Rpim': KEY_REFLNS_PDBX_RPIM_I_ALL, + 'CC½': KEY_REFLNS_PDBX_CC_HALF, +} + + +def _find_xia2_multiplex_overall_table(soup, panel_id='collapse_overall_All_data'): + """Find the "Overall" stats table for the combined ("All data") dataset. + + Returns a dict of row label -> [Overall, Low resolution, High resolution] cell text, + or None if the panel/table couldn't be found. + """ + panel = soup.find(id=panel_id) + if panel is None: + return None + table = panel.find('table') + if table is None: + return None + data = {} + for row in table.find_all('tr')[1:]: # first row is the column header + cells = row.find_all(['th', 'td']) + if not cells: + continue + data[cells[0].get_text(strip=True)] = [c.get_text(strip=True) for c in cells[1:]] + return data + + +def handle_xia2_multiplex(file): + """Parse the xia2.multiplex HTML report (not a plain-text log, unlike the other processing + programs), pulling stats from the "Overall" summary table of the combined "All data" dataset. + """ + reflns = {KEY_REFLNS_ENTRY_ID: 'UNNAMED', KEY_REFLNS_DIFFRN_ID: 1, KEY_REFLNS_PDBX_ORDINAL: 1} + shell = {KEY_REFLNS_DIFFRN_ID: (1, 1), KEY_REFLNS_PDBX_ORDINAL: (1, 2)} + + if file is None: + error('Log file not defined') + return reflns, shell + elif not Path(file).is_file(): + error('Log file ' + str(file) + ' not present') + return reflns, shell + + with open(file, 'rt', encoding='utf-8') as f: + soup = BeautifulSoup(f, 'html.parser') + + table_data = _find_xia2_multiplex_overall_table(soup) + if table_data is None: + warn('could not find "Overall" stats table in xia2.multiplex report ' + str(file)) + return reflns, shell + + resolution = table_data.get('Resolution (Å)') + if resolution is None or len(resolution) != 3: + warn('key ' + KEY_REFLNS_RESO_LOW + '/' + KEY_REFLNS_RESO_HIGH + ' not found for reflns') + else: + overall_low, overall_high = (v.strip() for v in resolution[0].split('-')) + low_shell_low, low_shell_high = (v.strip() for v in resolution[1].split('-')) + high_shell_low, high_shell_high = (v.strip() for v in resolution[2].split('-')) + reflns[KEY_REFLNS_RESO_LOW] = overall_low + reflns[KEY_REFLNS_RESO_HIGH] = overall_high + # first row is outer (high resolution) shell, second is inner (low resolution) shell + shell[_replace_shell_key(KEY_REFLNS_RESO_LOW)] = (high_shell_low, low_shell_low) + shell[_replace_shell_key(KEY_REFLNS_RESO_HIGH)] = (high_shell_high, low_shell_high) + + for label, key in _XIA2_MULTIPLEX_ROW_KEYS.items(): + values = table_data.get(label) + if values is None or len(values) != 3: + warn('key ' + key + ' not found for reflns') + continue + reflns[key] = values[0] + shell[_replace_shell_key(key)] = (values[2], values[1]) + + completeness = table_data.get('Completeness') + if completeness is None or len(completeness) != 3: + warn('key ' + KEY_REFLNS_POSSIBLE_OBS + ' not found for reflns') + else: + reflns[KEY_REFLNS_POSSIBLE_OBS] = completeness[0].rstrip('%') + shell[_replace_shell_key(KEY_REFLNS_POSSIBLE_OBS)] = ( + completeness[2].rstrip('%'), + completeness[1].rstrip('%'), + ) + + return reflns, shell + + def handle_file(file, type, doc: cif.Document, outputfile: str): if type == 'autoproc': reflns, shell = handle_autoproc(file) @@ -183,6 +276,8 @@ def handle_file(file, type, doc: cif.Document, outputfile: str): reflns, shell = handle_autoproc_staraniso(file) elif type == 'xia_3dii': reflns, shell = handle_xia_3dii(file) + elif type == 'xia2-multiplex': + reflns, shell = handle_xia2_multiplex(file) else: info('Unsupported type: ' + type) return None diff --git a/test-data/processing_stats/autoproc_aimless.log b/test-data/processing_stats/autoproc_aimless.log new file mode 100644 index 0000000..2885e2c --- /dev/null +++ b/test-data/processing_stats/autoproc_aimless.log @@ -0,0 +1,24 @@ +# Real excerpt from an autoPROC aimless.log (A71EV2A-x0836, lb32627-66), trimmed to the +# summary section that handle_autoproc()/handle_text_file() scans for. +Summary data for Project: lb32627v225 Crystal: xA71EV2Ax08361 Dataset: 1 + + Overall InnerShell OuterShell +Low resolution limit 46.38 46.38 1.41 +High resolution limit 1.34 4.24 1.34 + +Rmerge (within I+/I-) 0.162 0.054 2.339 +Rmerge (all I+ and I-) 0.173 0.057 2.543 +Rmeas (within I+/I-) 0.189 0.063 2.727 +Rmeas (all I+ & I-) 0.186 0.061 2.737 +Rpim (within I+/I-) 0.097 0.032 1.397 +Rpim (all I+ & I-) 0.068 0.022 1.005 +Rmerge in top intensity bin 0.053 - - +Total number of observations 210486 7398 27769 +Total number unique 28848 1019 3860 +Mean((I)/sd(I)) 6.5 23.5 0.7 +Mn(I) half-set correlation CC(1/2) 0.997 0.998 0.376 +Completeness 92.5 99.4 85.4 +Multiplicity 7.3 7.3 7.2 +Mean(Chi^2) 0.75 0.73 0.62 + +Anomalous completeness 91.5 99.7 84.0 diff --git a/test-data/processing_stats/staraniso_alldata-unique.table1 b/test-data/processing_stats/staraniso_alldata-unique.table1 new file mode 100644 index 0000000..cb5fc3a --- /dev/null +++ b/test-data/processing_stats/staraniso_alldata-unique.table1 @@ -0,0 +1,49 @@ + + Spacegroup name C2 + Unit cell parameters 73.059 61.473 32.723 90.000 92.442 90.000 + Wavelength 0.92203 A + + Diffraction limits & principal axes of ellipsoid fitted to diffraction cut-off surface: + 1.742 1.0000 0.0000 0.0000 _a_* - 0.019 _c_* + 1.536 0.0000 1.0000 0.0000 _b_* + 1.580 0.0000 0.0000 1.0000 _c_* + + Number of active ice-rings within this resolution range = 0 + + Criteria used in determination of diffraction limits: + ----------------------------------------------------- + local(I/sigI) >= 1.20 + + Per-reflection cut-off Operational Resolution + ----------------------------------------------------------------- + I/sigma(I) >= 2.0 : 1.9911 A for 9967 reflections + I/sigma(I) >= 1.0 : 1.9107 A for 11270 reflections + I/sigma(I) >= 0.0 : 1.8377 A for 12647 reflections + all : 1.7861 A for 13779 reflections + + Overall InnerShell OuterShell + --------------------------------------------------------------------------- + Low resolution limit 36.496 36.496 1.739 + High resolution limit 1.530 4.913 1.530 + + + Rmerge (all I+ & I-) 0.086 0.065 0.307 + Rmerge (within I+/I-) 0.082 0.062 0.315 + Rmeas (all I+ & I-) 0.095 0.071 0.406 + Rmeas (within I+/I-) 0.100 0.074 0.445 + Rpim (all I+ & I-) 0.040 0.029 0.262 + Rpim (within I+/I-) 0.055 0.040 0.314 + Total number of observations 49377 4384 1312 + Total number unique 13779 689 689 + Mean(I)/sd(I) 6.3 13.7 1.4 + Completeness (spherical) 63.0 100.0 9.9 + Completeness (ellipsoidal) 74.4 100.0 19.1 + Multiplicity 3.6 6.4 1.9 + CC(1/2) 0.996 0.989 0.818 + + + Anomalous completeness (spherical) 48.9 99.5 6.5 + Anomalous completeness (ellipsoidal) 57.7 99.5 12.5 + Anomalous multiplicity 2.2 3.4 1.2 + CC(ano) -0.372 -0.278 0.032 + |DANO|/sd(DANO) 0.464 0.464 0.438 diff --git a/test-data/processing_stats/xia2_multiplex.html b/test-data/processing_stats/xia2_multiplex.html new file mode 100644 index 0000000..3780025 --- /dev/null +++ b/test-data/processing_stats/xia2_multiplex.html @@ -0,0 +1,255 @@ + + +
+| + + + + | Overall | + + + +Low resolution | + + + +High resolution | + + +
|---|---|---|---|
| Resolution (Å) | + + + +87.11 - 1.89 | + + + +87.20 - 5.13 | + + + +1.92 - 1.89 | + + +
| Observations | + + + +358302 | + + + +18446 | + + + +12297 | + + +
| Unique reflections | + + + +30629 | + + + +1585 | + + + +1503 | + + +
| Multiplicity | + + + +11.7 | + + + +11.6 | + + + +8.2 | + + +
| Completeness | + + + +100.00% | + + + +100.00% | + + + +99.60% | + + +
| Mean I/σ(I) | + + + +14.6 | + + + +49.5 | + + + +1.8 | + + +
| Rmerge | + + + +0.282 | + + + +0.086 | + + + +3.395 | + + +
| Rmeas | + + + +0.295 | + + + +0.090 | + + + +3.629 | + + +
| Rpim | + + + +0.084 | + + + +0.026 | + + + +1.246 | + + +
| CC½ | + + + +0.998 | + + + +0.998 | + + + +0.351 | + + +
no stats here
") + + reflns, shell = sps.handle_xia2_multiplex(str(html_file)) + + assert reflns == {"entry_id": "UNNAMED", "pdbx_diffrn_id": 1, "pdbx_ordinal": 1} + assert shell == {"pdbx_diffrn_id": (1, 1), "pdbx_ordinal": (1, 2)} + + +# --------------------------------------------------------------------------- +# handle_file — dispatch by type, and CIF document assembly +# --------------------------------------------------------------------------- + + +def test_handle_file_autoproc_creates_reflns_loop(): + doc = sps.handle_file(str(AUTOPROC_LOG), "autoproc", None, None) + block = doc[0] + assert block.find_pair("_reflns.d_resolution_high") == ("_reflns.d_resolution_high", "1.34") + shell_item = block.find_loop_item("_reflns_shell.d_res_high") + assert list(shell_item.loop.values)[list(shell_item.loop.tags).index("_reflns_shell.d_res_high")] == "1.34" + + +def test_handle_file_autoproc_staraniso_creates_reflns_loop(): + doc = sps.handle_file(str(AUTOPROC_STARANISO_TABLE1), "autoproc_staraniso", None, None) + block = doc[0] + assert block.find_pair("_reflns.d_resolution_high") == ("_reflns.d_resolution_high", "1.530") + + +def test_handle_file_xia_3dii_creates_reflns_loop(): + doc = sps.handle_file(str(XIA_3DII_LOG), "xia_3dii", None, None) + block = doc[0] + assert block.find_pair("_reflns.d_resolution_high") == ("_reflns.d_resolution_high", "1.60") + + +def test_handle_file_xia2_multiplex_creates_reflns_loop(): + doc = sps.handle_file(str(XIA2_MULTIPLEX_HTML), "xia2-multiplex", None, None) + block = doc[0] + assert block.find_pair("_reflns.d_resolution_high") == ("_reflns.d_resolution_high", "1.89") + assert block.find_pair("_reflns.percent_possible_obs") == ("_reflns.percent_possible_obs", "100.00") + + +def test_handle_file_appends_to_existing_doc(): + """handle_file should add its new block to an existing doc rather than replacing it.""" + existing_doc = cif.Document() + existing_doc.add_new_block("existing") + + doc = sps.handle_file(str(XIA2_MULTIPLEX_HTML), "xia2-multiplex", existing_doc, None) + + assert len(doc) == 2 + assert doc[0].name == "existing" + + +def test_handle_file_unsupported_type_returns_none(): + assert sps.handle_file(str(AUTOPROC_LOG), "not-a-real-type", None, None) is None + + +def test_handle_file_writes_output_file(tmp_path): + out_file = tmp_path / "stats.cif" + sps.handle_file(str(XIA2_MULTIPLEX_HTML), "xia2-multiplex", None, str(out_file)) + assert out_file.is_file() + content = out_file.read_text() + assert "_reflns.d_resolution_high" in content