From b6c7e802d65efcc89b069e65bc4cffae1bba154c Mon Sep 17 00:00:00 2001 From: 3dot141592 Date: Tue, 21 Apr 2026 15:50:28 +0200 Subject: [PATCH 01/22] add geometry_operations.py and trimesh --- .../data_analysis/geometry_operations.py | 186 ++++++++++++++++++ requirements.txt | 2 + 2 files changed, 188 insertions(+) create mode 100644 backend/protzilla/data_analysis/geometry_operations.py diff --git a/backend/protzilla/data_analysis/geometry_operations.py b/backend/protzilla/data_analysis/geometry_operations.py new file mode 100644 index 000000000..379a77280 --- /dev/null +++ b/backend/protzilla/data_analysis/geometry_operations.py @@ -0,0 +1,186 @@ +import numpy as np +import pandas as pd +import trimesh +from trimesh import Trimesh +from trimesh.collision import CollisionManager + + +def _resolve_chain_column(cif_df: pd.DataFrame) -> str | None: + """ + Return the preferred chain identifier column if present in the CIF DataFrame. + """ + + if "_atom_site.label_asym_id" in cif_df.columns: + return "_atom_site.label_asym_id" + if "_atom_site.auth_asym_id" in cif_df.columns: + return "_atom_site.auth_asym_id" + return None + + +def extract_points_from_cif( + cif_df: pd.DataFrame, + residue_range: tuple[int, int], + chain_id: str | None = None, +) -> np.ndarray: + """ + Extract Cartesian atom coordinates in a residue range from a CIF DataFrame. + + :param cif_df: DataFrame containing mmCIF atom_site coordinates + :param residue_range: inclusive residue interval as [start, end] + :param chain_id: optional chain identifier. If given, the function will filter by + _atom_site.label_asym_id or _atom_site.auth_asym_id when available. + :return: array of n points und x, y, z coordinates. (n, 3) + :raises ValueError: if required columns are missing or too few points remain + """ + + required_columns = { + "_atom_site.label_seq_id", + "_atom_site.Cartn_x", + "_atom_site.Cartn_y", + "_atom_site.Cartn_z", + } + missing_columns = sorted(required_columns - set(cif_df.columns)) + if missing_columns: + raise ValueError( + f"CIF DataFrame is missing required columns for 3D extraction: {missing_columns}" + ) + + start, end = residue_range + if start > end: + raise ValueError( + f"Invalid residue range {residue_range}. Start must be smaller than or equal to end." + ) + + filtered_df = cif_df.copy() + + chain_column = _resolve_chain_column(filtered_df) + if chain_id is not None: + if chain_column is None: + raise ValueError( + "A chain_id was provided, but the CIF DataFrame has no chain identifier column." + ) + filtered_df = filtered_df[filtered_df[chain_column] == chain_id] + + residue_ids = filtered_df["_atom_site.label_seq_id"].astype(int) + filtered_df = filtered_df[(residue_ids >= start) & (residue_ids <= end)] + + points = ( + filtered_df[ + [ + "_atom_site.Cartn_x", + "_atom_site.Cartn_y", + "_atom_site.Cartn_z", + ] + ] + .astype(float) + .drop_duplicates() + .to_numpy() + ) + + if len(points) == 0: + raise ValueError("No atom coordinates found.") + + return points + + +def build_convex_hull(points: np.ndarray) -> Trimesh: + """ + Build a 3D convex hull mesh from a point cloud. + """ + + if points.ndim != 2 or points.shape[1] != 3: + raise ValueError( + f"Expected points with shape (n, 3), got array with shape {points.shape}." + ) + + if len(points) < 4: + raise ValueError( + "At least four distinct points are required to build a 3D convex hull." + ) + + return trimesh.convex.convex_hull(points, qhull_options="QJ") # Maybe QJ is stupid here? Ill have to look into it + + +def meshes_intersect( + mesh_a: Trimesh, mesh_b: Trimesh, distance_tolerance: float = 1e-9 +) -> bool: + """ + Determine whether two triangle meshes intersect or touch. + """ + + manager = CollisionManager() + manager.add_object("mesh_a", mesh_a) + return manager.min_distance_single(mesh_b) <= distance_tolerance + + +def meshes_distance(mesh_a: Trimesh, mesh_b: Trimesh) -> float: + """ + Calculate the minimum euclidean distance between two triangle meshes. + """ + + manager = CollisionManager() + manager.add_object("mesh_a", mesh_a) + distance = float(manager.min_distance_single(mesh_b)) + return max(distance, 0.0) + + +def bodies_intersect_from_cif( + cif_df: pd.DataFrame, + residue_range_a: tuple[int, int], + residue_range_b: tuple[int, int], + chain_id: str | None = None, +) -> dict: + """ + Build two convex bodies from CIF residue ranges and test whether they intersect. + + :param cif_df: DataFrame containing mmCIF atom_site coordinates + :param residue_range_a: inclusive residue range for the first body + :param residue_range_b: inclusive residue range for the second body + :param chain_id: optional chain identifier used for both bodies + :return: summary dictionary with hull sizes and the intersection result + """ + + points_a = extract_points_from_cif(cif_df, residue_range_a, chain_id=chain_id) + points_b = extract_points_from_cif(cif_df, residue_range_b, chain_id=chain_id) + + hull_a = build_convex_hull(points_a) + hull_b = build_convex_hull(points_b) + + return { + "intersects": meshes_intersect(hull_a, hull_b), + "n_atoms_a": len(points_a), + "n_atoms_b": len(points_b), + "n_hull_vertices_a": len(hull_a.vertices), + "n_hull_vertices_b": len(hull_b.vertices), + } + + +def bodies_distance_from_cif( + cif_df: pd.DataFrame, + residue_range_a: tuple[int, int], + residue_range_b: tuple[int, int], + chain_id: str | None = None, +) -> dict: + """ + Build two convex bodies from CIF residue ranges and calculate their distance. + + :param cif_df: DataFrame containing mmCIF atom_site coordinates + :param residue_range_a: inclusive residue range for the first body + :param residue_range_b: inclusive residue range for the second body + :param chain_id: optional chain identifier used for both bodies + :return: summary dictionary with hull sizes and the minimum distance + """ + + points_a = extract_points_from_cif(cif_df, residue_range_a, chain_id=chain_id) + points_b = extract_points_from_cif(cif_df, residue_range_b, chain_id=chain_id) + + hull_a = build_convex_hull(points_a) + hull_b = build_convex_hull(points_b) + + return { + "distance": meshes_distance(hull_a, hull_b), + "n_atoms_a": len(points_a), + "n_atoms_b": len(points_b), + "n_hull_vertices_a": len(hull_a.vertices), + "n_hull_vertices_b": len(hull_b.vertices), + } diff --git a/requirements.txt b/requirements.txt index b9119f4a9..992012ba4 100644 --- a/requirements.txt +++ b/requirements.txt @@ -32,6 +32,7 @@ pytest-django==4.11.1 pytest-order==1.3.0 PyYAML==6.0.3 py7zr!=1.1.2 +python-fcl==0.7.0.11 requests==2.32.5 restring==0.1.21 scikit-learn==1.7.2 @@ -40,6 +41,7 @@ seaborn==0.13.2 sphinx==8.2.3 sphinx-autoapi==3.6.1 statsmodels==0.14.5 +trimesh==4.11.5 tqdm==4.67.1 tzlocal==5.3.1 umap-learn==0.5.9.post2 From 1b3f70ad34d6694538314af9939369e6bc523cde Mon Sep 17 00:00:00 2001 From: 3dot141592 Date: Tue, 21 Apr 2026 16:24:01 +0200 Subject: [PATCH 02/22] run black --- backend/protzilla/data_analysis/geometry_operations.py | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/backend/protzilla/data_analysis/geometry_operations.py b/backend/protzilla/data_analysis/geometry_operations.py index 379a77280..8fa702e3d 100644 --- a/backend/protzilla/data_analysis/geometry_operations.py +++ b/backend/protzilla/data_analysis/geometry_operations.py @@ -98,7 +98,9 @@ def build_convex_hull(points: np.ndarray) -> Trimesh: "At least four distinct points are required to build a 3D convex hull." ) - return trimesh.convex.convex_hull(points, qhull_options="QJ") # Maybe QJ is stupid here? Ill have to look into it + return trimesh.convex.convex_hull( + points, qhull_options="QJ" + ) # Maybe QJ is stupid here? Ill have to look into it def meshes_intersect( From ea4f913b630ae4d0956413d8e05a8c48cde19bab Mon Sep 17 00:00:00 2001 From: 3dot141592 Date: Fri, 24 Apr 2026 14:09:01 +0200 Subject: [PATCH 03/22] LLM generated: Polyhedron visualizer --- .codex | 0 backend/main/views_helper.py | 8 ++ .../data_analysis/geometry_operations.py | 59 ++++++++++- backend/protzilla/run.py | 69 +++++++------ .../components/app/run-screen/run-screen.tsx | 10 +- .../core/shared/molstar-viewer/index.ts | 1 + .../molstar-viewer/molstar-trimesh-adapter.ts | 99 +++++++++++++++++++ .../molstar-viewer/molstar-viewer.props.tsx | 2 + .../shared/molstar-viewer/molstar-viewer.tsx | 15 ++- frontend/src/utils/protzilla-types.ts | 6 ++ 10 files changed, 229 insertions(+), 40 deletions(-) create mode 100644 .codex create mode 100644 frontend/src/components/core/shared/molstar-viewer/molstar-trimesh-adapter.ts diff --git a/.codex b/.codex new file mode 100644 index 000000000..e69de29bb diff --git a/backend/main/views_helper.py b/backend/main/views_helper.py index e7f6f7c0b..a3064bf53 100644 --- a/backend/main/views_helper.py +++ b/backend/main/views_helper.py @@ -6,6 +6,9 @@ from typing import Optional, List, Dict from backend.protzilla.constants.paths import SETTINGS_PATH +from backend.protzilla.data_analysis.geometry_operations import ( + convex_hull_polyhedron_from_cif, +) from backend.protzilla.disk_operator import YamlOperator from backend.protzilla.steps import Step from backend.protzilla.step_manager import StepManager @@ -209,6 +212,11 @@ def create_visualization( result = {"structureEntryId": structure_entry_id, "cifString": cif_string} + try: + result["polyhedron"] = convex_hull_polyhedron_from_cif(cif_df) + except (ValueError, TypeError, ModuleNotFoundError, ImportError): + pass + if crosslinking_df is not None: result["crosslinks"] = extract_relevant_crosslink_information(crosslinking_df) diff --git a/backend/protzilla/data_analysis/geometry_operations.py b/backend/protzilla/data_analysis/geometry_operations.py index 8fa702e3d..cb94de0e1 100644 --- a/backend/protzilla/data_analysis/geometry_operations.py +++ b/backend/protzilla/data_analysis/geometry_operations.py @@ -1,8 +1,12 @@ +from __future__ import annotations + import numpy as np import pandas as pd -import trimesh -from trimesh import Trimesh -from trimesh.collision import CollisionManager +from scipy.spatial import ConvexHull +from typing import TYPE_CHECKING + +if TYPE_CHECKING: + from trimesh import Trimesh def _resolve_chain_column(cif_df: pd.DataFrame) -> str | None: @@ -98,11 +102,56 @@ def build_convex_hull(points: np.ndarray) -> Trimesh: "At least four distinct points are required to build a 3D convex hull." ) + import trimesh + return trimesh.convex.convex_hull( points, qhull_options="QJ" ) # Maybe QJ is stupid here? Ill have to look into it +def convex_hull_polyhedron_from_cif(cif_df: pd.DataFrame) -> dict: + """ + Build a convex hull from all atom coordinates in a CIF DataFrame + and return it as triangle mesh lists for frontend visualization. + """ + + required_columns = [ + "_atom_site.Cartn_x", + "_atom_site.Cartn_y", + "_atom_site.Cartn_z", + ] + missing_columns = [column for column in required_columns if column not in cif_df.columns] + if missing_columns: + raise ValueError( + f"CIF DataFrame is missing required columns for convex hull creation: {missing_columns}" + ) + + coordinates = ( + cif_df[required_columns] + .apply(pd.to_numeric, errors="coerce") + .dropna() + .drop_duplicates() + .to_numpy() + ) + + try: + hull = build_convex_hull(coordinates) + return { + "vertices": hull.vertices.tolist(), + "faces": hull.faces.tolist(), + } + except ModuleNotFoundError: + hull = ConvexHull(coordinates) + used_vertex_indices = np.unique(hull.simplices) + remap = {old: new for new, old in enumerate(used_vertex_indices)} + faces = [[remap[index] for index in face] for face in hull.simplices.tolist()] + vertices = coordinates[used_vertex_indices] + return { + "vertices": vertices.tolist(), + "faces": faces, + } + + def meshes_intersect( mesh_a: Trimesh, mesh_b: Trimesh, distance_tolerance: float = 1e-9 ) -> bool: @@ -110,6 +159,8 @@ def meshes_intersect( Determine whether two triangle meshes intersect or touch. """ + from trimesh.collision import CollisionManager + manager = CollisionManager() manager.add_object("mesh_a", mesh_a) return manager.min_distance_single(mesh_b) <= distance_tolerance @@ -120,6 +171,8 @@ def meshes_distance(mesh_a: Trimesh, mesh_b: Trimesh) -> float: Calculate the minimum euclidean distance between two triangle meshes. """ + from trimesh.collision import CollisionManager + manager = CollisionManager() manager.add_object("mesh_a", mesh_a) distance = float(manager.min_distance_single(mesh_b)) diff --git a/backend/protzilla/run.py b/backend/protzilla/run.py index 3231f68ba..9cfc98d6b 100644 --- a/backend/protzilla/run.py +++ b/backend/protzilla/run.py @@ -63,38 +63,45 @@ def get_available_run_info() -> ( runs_favourited = [] all_tags = set() for run_name in get_available_run_names(): + try: + run_dir = os.path.join(paths.RUNS_PATH, run_name) + metadata_yaml_path = os.path.join(run_dir, "metadata.yaml") + if not os.path.isfile(metadata_yaml_path): + Run( + run_name + ) # initialize run to create metadata.yaml (creation date set to now) + yaml_operator = YamlOperator() + metadata = yaml_operator.read(Path(metadata_yaml_path)) + if not metadata: + metadata = {} + tags = metadata.get("tags", set()) + + run_name = { + "run_name": run_name, + "creation_date": metadata.get("creation_date", "date not available"), + "modification_date": metadata.get( + "modification_date", "date not available" + ), + "memory_mode": metadata.get("df_mode", "disk"), + "run_steps": metadata.get("steps", []), + "favourite_status": metadata.get("favourite", False), + "run_tags": list(tags), + } - run_dir = os.path.join(paths.RUNS_PATH, run_name) - metadata_yaml_path = os.path.join(run_dir, "metadata.yaml") - if not os.path.isfile(metadata_yaml_path): - Run( - run_name - ) # initialize run to create metadata.yaml (creation date set to now) - yaml_operator = YamlOperator() - metadata = yaml_operator.read(Path(metadata_yaml_path)) - if not metadata: - metadata = {} - tags = metadata.get("tags", set()) - - run_name = { - "run_name": run_name, - "creation_date": metadata.get("creation_date", "date not available"), - "modification_date": metadata.get( - "modification_date", "date not available" - ), - "memory_mode": metadata.get("df_mode", "disk"), - "run_steps": metadata.get("steps", []), - "favourite_status": metadata.get("favourite", False), - "run_tags": list(tags), - } - - if run_name["favourite_status"]: - runs_favourited.append(run_name) - else: - runs.append(run_name) - - for tag in tags: - all_tags.add(tag) + if run_name["favourite_status"]: + runs_favourited.append(run_name) + else: + runs.append(run_name) + + for tag in tags: + all_tags.add(tag) + except Exception as error: + logging.warning( + "Skipping invalid run folder '%s' while listing runs: %s", + run_name, + error, + ) + continue all_tags = list(all_tags) diff --git a/frontend/src/components/app/run-screen/run-screen.tsx b/frontend/src/components/app/run-screen/run-screen.tsx index cee1c0d28..5f150084f 100644 --- a/frontend/src/components/app/run-screen/run-screen.tsx +++ b/frontend/src/components/app/run-screen/run-screen.tsx @@ -30,7 +30,6 @@ import { Col } from "react-grid-system"; import { useLocation, useNavigate } from "react-router-dom"; import { styled } from "styled-components"; -import { CrosslinkerInformation } from "../../core/shared/molstar-viewer/crosslinker-processing"; import { H3 } from "../../core/shared/text"; const StyledNavbar = styled(Navbar)` @@ -162,13 +161,14 @@ export const RunScreen: React.FC = () => { structureEntryId: response.data.structureEntryId, cifString: response.data.cifString, crosslinks: response.data.crosslinks, + polyhedron: response.data.polyhedron, }), [], ); const visualizations = useCertainStepOutputs< StepOutputInfo, ApiResponse, - { structureEntryId: string; cifString: string; crosslinks?: CrosslinkerInformation[] } + Visualization >({ available_outputs: availableVisualizations, endpoint: "get_step_visualizations/", @@ -387,7 +387,11 @@ export const RunScreen: React.FC = () => { {visualizations.length > 0 ? ( visualizations.map((viz) => ( - + )) ) : ( diff --git a/frontend/src/components/core/shared/molstar-viewer/index.ts b/frontend/src/components/core/shared/molstar-viewer/index.ts index 45fc7b6db..005c3ac68 100644 --- a/frontend/src/components/core/shared/molstar-viewer/index.ts +++ b/frontend/src/components/core/shared/molstar-viewer/index.ts @@ -1,2 +1,3 @@ export { default as MolstarViewer } from "./molstar-viewer"; export * from "./molstar-viewer.props"; +export * from "./molstar-trimesh-adapter"; diff --git a/frontend/src/components/core/shared/molstar-viewer/molstar-trimesh-adapter.ts b/frontend/src/components/core/shared/molstar-viewer/molstar-trimesh-adapter.ts new file mode 100644 index 000000000..de5baf5f2 --- /dev/null +++ b/frontend/src/components/core/shared/molstar-viewer/molstar-trimesh-adapter.ts @@ -0,0 +1,99 @@ +import { Mesh } from "molstar/lib/mol-geo/geometry/mesh/mesh"; +import { MeshBuilder } from "molstar/lib/mol-geo/geometry/mesh/mesh-builder"; +import { Vec3 } from "molstar/lib/mol-math/linear-algebra"; +import { Shape } from "molstar/lib/mol-model/shape"; +import { PluginStateObject } from "molstar/lib/mol-plugin-state/objects"; +import { StateTransforms } from "molstar/lib/mol-plugin-state/transforms"; +import { PluginUIContext } from "molstar/lib/mol-plugin-ui/context"; +import { StateTransformer } from "molstar/lib/mol-state"; +import { Task } from "molstar/lib/mol-task"; +import { Color } from "molstar/lib/mol-util/color"; +import { ParamDefinition as PD } from "molstar/lib/mol-util/param-definition"; + +export interface TrimeshLike { + vertices: number[][]; + faces: number[][]; +} + +interface AddPolyhedronOptions { + color?: number; + alpha?: number; + label?: string; +} + +const ProtzillaTransforms = StateTransformer.builderFactory("protzilla"); + +const TrimeshShape = ProtzillaTransforms({ + name: "trimesh-shape", + display: { name: "Trimesh Shape" }, + from: PluginStateObject.Root, + to: PluginStateObject.Shape.Provider, + params: { + mesh: PD.Value({ vertices: [], faces: [] }, { isHidden: true }), + color: PD.Color(Color(0xff8800)), + label: PD.Text("Polyhedron"), + }, +})({ + canAutoUpdate: () => true, + apply({ params }) { + return Task.create("Create Trimesh Shape", () => { + return new PluginStateObject.Shape.Provider( + { + label: params.label, + data: params.mesh, + params: Mesh.Params, + geometryUtils: Mesh.Utils, + getShape: (_ctx, data) => { + const builder = MeshBuilder.createState(data.vertices.length, data.faces.length); + builder.currentGroup = 0; + + for (const face of data.faces) { + if (face.length !== 3) continue; + const a = data.vertices[face[0]]; + const b = data.vertices[face[1]]; + const c = data.vertices[face[2]]; + if (a === undefined || b === undefined || c === undefined) continue; + + MeshBuilder.addTriangle( + builder, + Vec3.create(a[0], a[1], a[2]), + Vec3.create(b[0], b[1], b[2]), + Vec3.create(c[0], c[1], c[2]), + ); + } + + const mesh = MeshBuilder.getMesh(builder); + return Shape.create( + params.label, + data, + mesh, + () => params.color, + () => 1, + () => params.label, + ); + }, + }, + { label: params.label }, + ); + }); + }, +}); + +export async function addTrimeshPolyhedron( + plugin: PluginUIContext, + mesh: TrimeshLike, + options: AddPolyhedronOptions = {}, +) { + if (mesh.vertices.length === 0 || mesh.faces.length === 0) return; + + const label = options.label ?? "Polyhedron"; + const color = Color(options.color ?? 0xff8800); + const alpha = options.alpha ?? 0.35; + + await plugin + .build() + .toRoot() + .apply(TrimeshShape, { mesh, color, label }) + .apply(StateTransforms.Representation.ShapeRepresentation3D, { alpha }) + .commit(); +} diff --git a/frontend/src/components/core/shared/molstar-viewer/molstar-viewer.props.tsx b/frontend/src/components/core/shared/molstar-viewer/molstar-viewer.props.tsx index abf936f4d..083c360aa 100644 --- a/frontend/src/components/core/shared/molstar-viewer/molstar-viewer.props.tsx +++ b/frontend/src/components/core/shared/molstar-viewer/molstar-viewer.props.tsx @@ -1,6 +1,8 @@ import { CrosslinkerInformation } from "./crosslinker-processing"; +import { TrimeshLike } from "./molstar-trimesh-adapter"; export interface MolstarViewerProps { cifText: string; crosslinks?: CrosslinkerInformation[]; + polyhedron?: TrimeshLike; } diff --git a/frontend/src/components/core/shared/molstar-viewer/molstar-viewer.tsx b/frontend/src/components/core/shared/molstar-viewer/molstar-viewer.tsx index 34ee45713..7c14b65b4 100644 --- a/frontend/src/components/core/shared/molstar-viewer/molstar-viewer.tsx +++ b/frontend/src/components/core/shared/molstar-viewer/molstar-viewer.tsx @@ -3,14 +3,15 @@ import { SectionTitle } from "@protzilla/core"; import { createPluginUI } from "molstar/lib/mol-plugin-ui"; import { PluginUIContext } from "molstar/lib/mol-plugin-ui/context"; import { renderReact18 } from "molstar/lib/mol-plugin-ui/react18"; +import "molstar/lib/mol-plugin-ui/skin/base/base.scss"; import React, { useEffect, useRef, useState } from "react"; +import { addTrimeshPolyhedron } from "./molstar-trimesh-adapter"; import { MolstarViewerProps } from "./molstar-viewer.props"; import { addCrosslinks, handleError } from "./molstar-viewer.service"; import { CanvasWrapper, Container } from "./styles"; -import "molstar/lib/mol-plugin-ui/skin/base/base.scss"; -const MolstarViewer: React.FC = ({ cifText, crosslinks }) => { +const MolstarViewer: React.FC = ({ cifText, crosslinks, polyhedron }) => { const containerRef = useRef(null); const [isLoading, setIsLoading] = useState(true); const notify = useNotification(); @@ -48,6 +49,14 @@ const MolstarViewer: React.FC = ({ cifText, crosslinks }) => await addCrosslinks(plugin, cifText, crosslinks); } + if (polyhedron !== undefined) { + await addTrimeshPolyhedron(plugin, polyhedron, { + color: 0x8a2be2, + alpha: 0.5, + label: "Protein Convex Hull", + }); + } + setIsLoading(false); } catch (error: unknown) { handleError(error, "MolstarViewer Error:", notify); @@ -66,7 +75,7 @@ const MolstarViewer: React.FC = ({ cifText, crosslinks }) => } } }; - }, [cifText, crosslinks, notify]); + }, [cifText, crosslinks, notify, polyhedron]); return ( diff --git a/frontend/src/utils/protzilla-types.ts b/frontend/src/utils/protzilla-types.ts index 4317559bd..cb308be27 100644 --- a/frontend/src/utils/protzilla-types.ts +++ b/frontend/src/utils/protzilla-types.ts @@ -35,10 +35,16 @@ export interface Download { data: Record; } +export interface PolyhedronData { + vertices: number[][]; + faces: number[][]; +} + export interface Visualization { structureEntryId: string; cifString: string; crosslinks?: CrosslinkerInformation[]; + polyhedron?: PolyhedronData; } // We assume these are the only data types we receive for tables From ff9244163bdecc94dd3b1eb366daed2b58c6806d Mon Sep 17 00:00:00 2001 From: 3dot141592 Date: Wed, 6 May 2026 16:50:23 +0200 Subject: [PATCH 04/22] Render amino acid spheres --- .codex | 0 backend/main/views_helper.py | 13 +- .../data_analysis/amino_acid_spheres.py | 53 ++++++ .../data_analysis/geometry_operations.py | 153 +++++------------- .../components/app/run-screen/run-screen.tsx | 4 +- .../molstar-viewer/molstar-trimesh-adapter.ts | 90 ++++++----- .../molstar-viewer/molstar-viewer.props.tsx | 11 +- .../shared/molstar-viewer/molstar-viewer.tsx | 20 +-- frontend/src/utils/protzilla-types.ts | 11 +- 9 files changed, 176 insertions(+), 179 deletions(-) delete mode 100644 .codex create mode 100644 backend/protzilla/data_analysis/amino_acid_spheres.py diff --git a/.codex b/.codex deleted file mode 100644 index e69de29bb..000000000 diff --git a/backend/main/views_helper.py b/backend/main/views_helper.py index a3064bf53..29e5f6e9d 100644 --- a/backend/main/views_helper.py +++ b/backend/main/views_helper.py @@ -6,9 +6,7 @@ from typing import Optional, List, Dict from backend.protzilla.constants.paths import SETTINGS_PATH -from backend.protzilla.data_analysis.geometry_operations import ( - convex_hull_polyhedron_from_cif, -) +from backend.protzilla.data_analysis.amino_acid_spheres import calculate_amino_acid_spheres from backend.protzilla.disk_operator import YamlOperator from backend.protzilla.steps import Step from backend.protzilla.step_manager import StepManager @@ -211,12 +209,9 @@ def create_visualization( cif_string = "" result = {"structureEntryId": structure_entry_id, "cifString": cif_string} - - try: - result["polyhedron"] = convex_hull_polyhedron_from_cif(cif_df) - except (ValueError, TypeError, ModuleNotFoundError, ImportError): - pass - + + result["trimeshMeshes"] = calculate_amino_acid_spheres(cif_df) + if crosslinking_df is not None: result["crosslinks"] = extract_relevant_crosslink_information(crosslinking_df) diff --git a/backend/protzilla/data_analysis/amino_acid_spheres.py b/backend/protzilla/data_analysis/amino_acid_spheres.py new file mode 100644 index 000000000..c9cd57ace --- /dev/null +++ b/backend/protzilla/data_analysis/amino_acid_spheres.py @@ -0,0 +1,53 @@ +from __future__ import annotations + +import numpy as np +import pandas as pd +import trimesh + +from backend.protzilla.data_analysis.geometry_operations import ( + calculate_centroid, + extract_points_from_cif, + find_farthest_point, + mesh_to_polyhedron, +) + + +def calculate_amino_acid_spheres( + cif_df: pd.DataFrame, + chain_id: str | None = None, + color: int = 0xFF8C00, + alpha: float = 0.25, + subdivisions: int = 1, +) -> list[dict]: + residue_positions = ( + pd.to_numeric(cif_df["_atom_site.label_seq_id"], errors="coerce") + .dropna() + .astype(int) + .drop_duplicates() + .sort_values() + ) + + spheres = [] + for residue_position in residue_positions: # TODO: This can certainly be made more efficient + residue_points = extract_points_from_cif( + cif_df, + residue_range=(residue_position, residue_position), + chain_id=chain_id, + ) + center = calculate_centroid(residue_points) + furthest_point = find_farthest_point(residue_points, center) + radius = float(np.linalg.norm(furthest_point - center)) + + sphere = trimesh.creation.icosphere(subdivisions=subdivisions, radius=radius) + sphere.apply_translation(center) + + spheres.append( + { + "label": f"Residue {residue_position} sphere", + "mesh": mesh_to_polyhedron(sphere), + "color": color, + "alpha": alpha, + } + ) + + return spheres diff --git a/backend/protzilla/data_analysis/geometry_operations.py b/backend/protzilla/data_analysis/geometry_operations.py index cb94de0e1..7b22e29cd 100644 --- a/backend/protzilla/data_analysis/geometry_operations.py +++ b/backend/protzilla/data_analysis/geometry_operations.py @@ -2,11 +2,13 @@ import numpy as np import pandas as pd -from scipy.spatial import ConvexHull -from typing import TYPE_CHECKING +import trimesh -if TYPE_CHECKING: - from trimesh import Trimesh +COORDINATE_COLUMNS = [ + "_atom_site.Cartn_x", + "_atom_site.Cartn_y", + "_atom_site.Cartn_z", +] def _resolve_chain_column(cif_df: pd.DataFrame) -> str | None: @@ -87,7 +89,7 @@ def extract_points_from_cif( return points -def build_convex_hull(points: np.ndarray) -> Trimesh: +def build_convex_hull(points: np.ndarray) -> trimesh.Trimesh: """ Build a 3D convex hull mesh from a point cloud. """ @@ -102,61 +104,50 @@ def build_convex_hull(points: np.ndarray) -> Trimesh: "At least four distinct points are required to build a 3D convex hull." ) - import trimesh + return trimesh.convex.convex_hull(points, qhull_options="QJ") - return trimesh.convex.convex_hull( - points, qhull_options="QJ" - ) # Maybe QJ is stupid here? Ill have to look into it +def mesh_to_polyhedron(mesh: trimesh.Trimesh) -> dict: + return { + "vertices": mesh.vertices.tolist(), + "faces": mesh.faces.tolist(), + } -def convex_hull_polyhedron_from_cif(cif_df: pd.DataFrame) -> dict: - """ - Build a convex hull from all atom coordinates in a CIF DataFrame - and return it as triangle mesh lists for frontend visualization. - """ - required_columns = [ - "_atom_site.Cartn_x", - "_atom_site.Cartn_y", - "_atom_site.Cartn_z", - ] - missing_columns = [column for column in required_columns if column not in cif_df.columns] - if missing_columns: +def point_cloud_to_polyhedron(points: np.ndarray) -> dict: + return mesh_to_polyhedron(build_convex_hull(points)) + + +def calculate_centroid(points: np.ndarray) -> np.ndarray: + if points.ndim != 2 or points.shape[1] != 3: + raise ValueError(f"Expected points with shape (n, 3), got {points.shape}.") + if len(points) == 0: + raise ValueError("At least one point is required to calculate a centroid.") + + return points.mean(axis=0) + + +def find_farthest_point(points: np.ndarray, reference_point: np.ndarray) -> np.ndarray: + if points.ndim != 2 or points.shape[1] != 3: + raise ValueError(f"Expected points with shape (n, 3), got {points.shape}.") + if len(points) == 0: + raise ValueError("At least one point is required to calculate a maximum distance.") + if reference_point.shape != (3,): raise ValueError( - f"CIF DataFrame is missing required columns for convex hull creation: {missing_columns}" + f"Expected reference_point with shape (3,), got {reference_point.shape}." ) - coordinates = ( - cif_df[required_columns] - .apply(pd.to_numeric, errors="coerce") - .dropna() - .drop_duplicates() - .to_numpy() - ) - - try: - hull = build_convex_hull(coordinates) - return { - "vertices": hull.vertices.tolist(), - "faces": hull.faces.tolist(), - } - except ModuleNotFoundError: - hull = ConvexHull(coordinates) - used_vertex_indices = np.unique(hull.simplices) - remap = {old: new for new, old in enumerate(used_vertex_indices)} - faces = [[remap[index] for index in face] for face in hull.simplices.tolist()] - vertices = coordinates[used_vertex_indices] - return { - "vertices": vertices.tolist(), - "faces": faces, - } + distances = np.linalg.norm(points - reference_point, axis=1) + return points[np.argmax(distances)] def meshes_intersect( - mesh_a: Trimesh, mesh_b: Trimesh, distance_tolerance: float = 1e-9 + mesh_a: trimesh.Trimesh, + mesh_b: trimesh.Trimesh, + distance_tolerance: float = 1e-9, ) -> bool: """ - Determine whether two triangle meshes intersect or touch. + Determine whether two meshes intersect or touch. """ from trimesh.collision import CollisionManager @@ -166,9 +157,9 @@ def meshes_intersect( return manager.min_distance_single(mesh_b) <= distance_tolerance -def meshes_distance(mesh_a: Trimesh, mesh_b: Trimesh) -> float: +def meshes_distance(mesh_a: trimesh.Trimesh, mesh_b: trimesh.Trimesh) -> float: """ - Calculate the minimum euclidean distance between two triangle meshes. + Calculate the minimum euclidean distance between two meshes. """ from trimesh.collision import CollisionManager @@ -177,65 +168,3 @@ def meshes_distance(mesh_a: Trimesh, mesh_b: Trimesh) -> float: manager.add_object("mesh_a", mesh_a) distance = float(manager.min_distance_single(mesh_b)) return max(distance, 0.0) - - -def bodies_intersect_from_cif( - cif_df: pd.DataFrame, - residue_range_a: tuple[int, int], - residue_range_b: tuple[int, int], - chain_id: str | None = None, -) -> dict: - """ - Build two convex bodies from CIF residue ranges and test whether they intersect. - - :param cif_df: DataFrame containing mmCIF atom_site coordinates - :param residue_range_a: inclusive residue range for the first body - :param residue_range_b: inclusive residue range for the second body - :param chain_id: optional chain identifier used for both bodies - :return: summary dictionary with hull sizes and the intersection result - """ - - points_a = extract_points_from_cif(cif_df, residue_range_a, chain_id=chain_id) - points_b = extract_points_from_cif(cif_df, residue_range_b, chain_id=chain_id) - - hull_a = build_convex_hull(points_a) - hull_b = build_convex_hull(points_b) - - return { - "intersects": meshes_intersect(hull_a, hull_b), - "n_atoms_a": len(points_a), - "n_atoms_b": len(points_b), - "n_hull_vertices_a": len(hull_a.vertices), - "n_hull_vertices_b": len(hull_b.vertices), - } - - -def bodies_distance_from_cif( - cif_df: pd.DataFrame, - residue_range_a: tuple[int, int], - residue_range_b: tuple[int, int], - chain_id: str | None = None, -) -> dict: - """ - Build two convex bodies from CIF residue ranges and calculate their distance. - - :param cif_df: DataFrame containing mmCIF atom_site coordinates - :param residue_range_a: inclusive residue range for the first body - :param residue_range_b: inclusive residue range for the second body - :param chain_id: optional chain identifier used for both bodies - :return: summary dictionary with hull sizes and the minimum distance - """ - - points_a = extract_points_from_cif(cif_df, residue_range_a, chain_id=chain_id) - points_b = extract_points_from_cif(cif_df, residue_range_b, chain_id=chain_id) - - hull_a = build_convex_hull(points_a) - hull_b = build_convex_hull(points_b) - - return { - "distance": meshes_distance(hull_a, hull_b), - "n_atoms_a": len(points_a), - "n_atoms_b": len(points_b), - "n_hull_vertices_a": len(hull_a.vertices), - "n_hull_vertices_b": len(hull_b.vertices), - } diff --git a/frontend/src/components/app/run-screen/run-screen.tsx b/frontend/src/components/app/run-screen/run-screen.tsx index 5f150084f..1b9ce1668 100644 --- a/frontend/src/components/app/run-screen/run-screen.tsx +++ b/frontend/src/components/app/run-screen/run-screen.tsx @@ -161,7 +161,7 @@ export const RunScreen: React.FC = () => { structureEntryId: response.data.structureEntryId, cifString: response.data.cifString, crosslinks: response.data.crosslinks, - polyhedron: response.data.polyhedron, + trimeshMeshes: response.data.trimeshMeshes, }), [], ); @@ -390,7 +390,7 @@ export const RunScreen: React.FC = () => { )) diff --git a/frontend/src/components/core/shared/molstar-viewer/molstar-trimesh-adapter.ts b/frontend/src/components/core/shared/molstar-viewer/molstar-trimesh-adapter.ts index de5baf5f2..859e060ca 100644 --- a/frontend/src/components/core/shared/molstar-viewer/molstar-trimesh-adapter.ts +++ b/frontend/src/components/core/shared/molstar-viewer/molstar-trimesh-adapter.ts @@ -10,12 +10,12 @@ import { Task } from "molstar/lib/mol-task"; import { Color } from "molstar/lib/mol-util/color"; import { ParamDefinition as PD } from "molstar/lib/mol-util/param-definition"; -export interface TrimeshLike { +export interface TrimeshMesh { vertices: number[][]; faces: number[][]; } -interface AddPolyhedronOptions { +export interface AddTrimeshMeshOptions { color?: number; alpha?: number; label?: string; @@ -23,77 +23,81 @@ interface AddPolyhedronOptions { const ProtzillaTransforms = StateTransformer.builderFactory("protzilla"); -const TrimeshShape = ProtzillaTransforms({ +const TrimeshMeshShape = ProtzillaTransforms({ name: "trimesh-shape", display: { name: "Trimesh Shape" }, from: PluginStateObject.Root, to: PluginStateObject.Shape.Provider, params: { - mesh: PD.Value({ vertices: [], faces: [] }, { isHidden: true }), + mesh: PD.Value({ vertices: [], faces: [] }, { isHidden: true }), color: PD.Color(Color(0xff8800)), - label: PD.Text("Polyhedron"), + label: PD.Text("Trimesh Mesh"), }, })({ canAutoUpdate: () => true, apply({ params }) { return Task.create("Create Trimesh Shape", () => { - return new PluginStateObject.Shape.Provider( - { - label: params.label, - data: params.mesh, - params: Mesh.Params, - geometryUtils: Mesh.Utils, - getShape: (_ctx, data) => { - const builder = MeshBuilder.createState(data.vertices.length, data.faces.length); - builder.currentGroup = 0; + return Promise.resolve( + new PluginStateObject.Shape.Provider( + { + label: params.label, + data: params.mesh, + params: Mesh.Params, + geometryUtils: Mesh.Utils, + getShape: (_ctx, data) => { + const builder = MeshBuilder.createState(data.vertices.length, data.faces.length); + builder.currentGroup = 0; - for (const face of data.faces) { - if (face.length !== 3) continue; - const a = data.vertices[face[0]]; - const b = data.vertices[face[1]]; - const c = data.vertices[face[2]]; - if (a === undefined || b === undefined || c === undefined) continue; + for (const face of data.faces) { + if (face.length !== 3) continue; + const a = data.vertices[face[0]]; + const b = data.vertices[face[1]]; + const c = data.vertices[face[2]]; + if (a === undefined || b === undefined || c === undefined) continue; - MeshBuilder.addTriangle( - builder, - Vec3.create(a[0], a[1], a[2]), - Vec3.create(b[0], b[1], b[2]), - Vec3.create(c[0], c[1], c[2]), - ); - } + MeshBuilder.addTriangle( + builder, + Vec3.create(a[0], a[1], a[2]), + Vec3.create(b[0], b[1], b[2]), + Vec3.create(c[0], c[1], c[2]), + ); + } - const mesh = MeshBuilder.getMesh(builder); - return Shape.create( - params.label, - data, - mesh, - () => params.color, - () => 1, - () => params.label, - ); + const mesh = MeshBuilder.getMesh(builder); + return Shape.create( + params.label, + data, + mesh, + () => params.color, + () => 1, + () => params.label, + ); + }, }, - }, - { label: params.label }, + { label: params.label }, + ), ); }); }, }); -export async function addTrimeshPolyhedron( +export async function addTrimeshMesh( plugin: PluginUIContext, - mesh: TrimeshLike, - options: AddPolyhedronOptions = {}, + mesh: TrimeshMesh, + options: AddTrimeshMeshOptions = {}, ) { if (mesh.vertices.length === 0 || mesh.faces.length === 0) return; - const label = options.label ?? "Polyhedron"; + const label = options.label ?? "Trimesh Mesh"; const color = Color(options.color ?? 0xff8800); const alpha = options.alpha ?? 0.35; await plugin .build() .toRoot() - .apply(TrimeshShape, { mesh, color, label }) + .apply(TrimeshMeshShape, { mesh, color, label }) .apply(StateTransforms.Representation.ShapeRepresentation3D, { alpha }) .commit(); } + +export const addTrimeshPolyhedron = addTrimeshMesh; diff --git a/frontend/src/components/core/shared/molstar-viewer/molstar-viewer.props.tsx b/frontend/src/components/core/shared/molstar-viewer/molstar-viewer.props.tsx index 083c360aa..c19f0391c 100644 --- a/frontend/src/components/core/shared/molstar-viewer/molstar-viewer.props.tsx +++ b/frontend/src/components/core/shared/molstar-viewer/molstar-viewer.props.tsx @@ -1,8 +1,15 @@ import { CrosslinkerInformation } from "./crosslinker-processing"; -import { TrimeshLike } from "./molstar-trimesh-adapter"; +import { TrimeshMesh } from "./molstar-trimesh-adapter"; + +export interface TrimeshShape { + label: string; + mesh: TrimeshMesh; + color?: number; + alpha?: number; +} export interface MolstarViewerProps { cifText: string; crosslinks?: CrosslinkerInformation[]; - polyhedron?: TrimeshLike; + trimeshMeshes?: TrimeshShape[]; } diff --git a/frontend/src/components/core/shared/molstar-viewer/molstar-viewer.tsx b/frontend/src/components/core/shared/molstar-viewer/molstar-viewer.tsx index 7c14b65b4..2f5ee65ec 100644 --- a/frontend/src/components/core/shared/molstar-viewer/molstar-viewer.tsx +++ b/frontend/src/components/core/shared/molstar-viewer/molstar-viewer.tsx @@ -6,12 +6,12 @@ import { renderReact18 } from "molstar/lib/mol-plugin-ui/react18"; import "molstar/lib/mol-plugin-ui/skin/base/base.scss"; import React, { useEffect, useRef, useState } from "react"; -import { addTrimeshPolyhedron } from "./molstar-trimesh-adapter"; +import { addTrimeshMesh } from "./molstar-trimesh-adapter"; import { MolstarViewerProps } from "./molstar-viewer.props"; import { addCrosslinks, handleError } from "./molstar-viewer.service"; import { CanvasWrapper, Container } from "./styles"; -const MolstarViewer: React.FC = ({ cifText, crosslinks, polyhedron }) => { +const MolstarViewer: React.FC = ({ cifText, crosslinks, trimeshMeshes }) => { const containerRef = useRef(null); const [isLoading, setIsLoading] = useState(true); const notify = useNotification(); @@ -49,12 +49,14 @@ const MolstarViewer: React.FC = ({ cifText, crosslinks, poly await addCrosslinks(plugin, cifText, crosslinks); } - if (polyhedron !== undefined) { - await addTrimeshPolyhedron(plugin, polyhedron, { - color: 0x8a2be2, - alpha: 0.5, - label: "Protein Convex Hull", - }); + if (trimeshMeshes !== undefined) { + for (const trimeshMesh of trimeshMeshes) { + await addTrimeshMesh(plugin, trimeshMesh.mesh, { + color: trimeshMesh.color, + alpha: trimeshMesh.alpha, + label: trimeshMesh.label, + }); + } } setIsLoading(false); @@ -75,7 +77,7 @@ const MolstarViewer: React.FC = ({ cifText, crosslinks, poly } } }; - }, [cifText, crosslinks, notify, polyhedron]); + }, [cifText, crosslinks, notify, trimeshMeshes]); return ( diff --git a/frontend/src/utils/protzilla-types.ts b/frontend/src/utils/protzilla-types.ts index cb308be27..c6d48ce5b 100644 --- a/frontend/src/utils/protzilla-types.ts +++ b/frontend/src/utils/protzilla-types.ts @@ -35,16 +35,23 @@ export interface Download { data: Record; } -export interface PolyhedronData { +export interface TrimeshMesh { vertices: number[][]; faces: number[][]; } +export interface TrimeshShape { + label: string; + mesh: TrimeshMesh; + color?: number; + alpha?: number; +} + export interface Visualization { structureEntryId: string; cifString: string; crosslinks?: CrosslinkerInformation[]; - polyhedron?: PolyhedronData; + trimeshMeshes?: TrimeshShape[]; } // We assume these are the only data types we receive for tables From ccb66daef6f13059658e1026e5cfb6a6f1f13655 Mon Sep 17 00:00:00 2001 From: Tarek Massini Date: Sun, 10 May 2026 20:55:46 +0200 Subject: [PATCH 05/22] feat: introduce parsing of _chem_comp table in cif-files --- backend/protzilla/constants/cif_columns.py | 63 +++++ .../alphafold_protein_structure_load.py | 69 +++-- .../test_alphafold_protein_structure_load.py | 260 ++++++++++++------ 3 files changed, 285 insertions(+), 107 deletions(-) create mode 100644 backend/protzilla/constants/cif_columns.py diff --git a/backend/protzilla/constants/cif_columns.py b/backend/protzilla/constants/cif_columns.py new file mode 100644 index 000000000..5715f2f32 --- /dev/null +++ b/backend/protzilla/constants/cif_columns.py @@ -0,0 +1,63 @@ +from enum import StrEnum + + +ATOM_SITE_PREFIX = "_atom_site." + + +class ATOM_SITE_COLUMNS(StrEnum): + """ + Enum containing all column names that should be present in + the _atom_site. table for mmCIF files from PDB or AFDB + """ + + ID = f"{ATOM_SITE_PREFIX}id" + TYPE_SYMBOL = f"{ATOM_SITE_PREFIX}type_symbol" + LABEL_ATOM_ID = f"{ATOM_SITE_PREFIX}label_atom_id" + LABEL_ALT_ID = f"{ATOM_SITE_PREFIX}label_alt_id" + LABEL_COMP_ID = f"{ATOM_SITE_PREFIX}label_comp_id" + LABEL_ASYM_ID = f"{ATOM_SITE_PREFIX}label_asym_id" + LABEL_ENTITY_ID = f"{ATOM_SITE_PREFIX}label_entity_id" + LABEL_SEQ_ID = f"{ATOM_SITE_PREFIX}label_seq_id" + PDBX_PDB_INS_CODE = f"{ATOM_SITE_PREFIX}pdbx_PDB_ins_code" + CARTN_X = f"{ATOM_SITE_PREFIX}Cartn_x" + CARTN_Y = f"{ATOM_SITE_PREFIX}Cartn_y" + CARTN_Z = f"{ATOM_SITE_PREFIX}Cartn_z" + OCCUPANCY = f"{ATOM_SITE_PREFIX}occupancy" + B_ISO_OR_EQUIV = f"{ATOM_SITE_PREFIX}B_iso_or_equiv" + PDBX_FORMAL_CHARGE = f"{ATOM_SITE_PREFIX}pdbx_formal_charge" + AUTH_SEQ_ID = f"{ATOM_SITE_PREFIX}auth_seq_id" + AUTH_COMP_ID = f"{ATOM_SITE_PREFIX}auth_comp_id" + AUTH_ASYM_ID = f"{ATOM_SITE_PREFIX}auth_asym_id" + AUTH_ATOM_ID = f"{ATOM_SITE_PREFIX}auth_atom_id" + PDBX_PDB_MODEL_NUM = f"{ATOM_SITE_PREFIX}pdbx_PDB_model_num" + + +ATOM_SITE_LABEL_COMP_ID = ATOM_SITE_COLUMNS.LABEL_COMP_ID + +ATOM_SITE_COLUMNS_NUMERIC = [ + ATOM_SITE_COLUMNS.ID, + ATOM_SITE_COLUMNS.LABEL_SEQ_ID, + ATOM_SITE_COLUMNS.CARTN_X, + ATOM_SITE_COLUMNS.CARTN_Y, + ATOM_SITE_COLUMNS.CARTN_Z, + ATOM_SITE_COLUMNS.OCCUPANCY, + ATOM_SITE_COLUMNS.B_ISO_OR_EQUIV, + ATOM_SITE_COLUMNS.AUTH_SEQ_ID, +] + +CHEM_COMP_PREFIX = "_chem_comp." + + +class CHEM_COMP_COLUMNS(StrEnum): + """ + Enum containing all column names that should be present in + the _chem_comp. table for mmCIF files from PDB or AFDB + """ + + ID = f"{CHEM_COMP_PREFIX}id" + TYPE = f"{CHEM_COMP_PREFIX}type" + MON_NSTD_FLAG = f"{CHEM_COMP_PREFIX}mon_nstd_flag" + NAME = f"{CHEM_COMP_PREFIX}name" + PDBX_SYNONYMS = f"{CHEM_COMP_PREFIX}pdbx_synonyms" + FORMULA = f"{CHEM_COMP_PREFIX}formula" + FORMULA_WEIGHT = f"{CHEM_COMP_PREFIX}formula_weight" diff --git a/backend/protzilla/importing/alphafold_protein_structure_load.py b/backend/protzilla/importing/alphafold_protein_structure_load.py index 7519c1d9b..49d9e10e6 100644 --- a/backend/protzilla/importing/alphafold_protein_structure_load.py +++ b/backend/protzilla/importing/alphafold_protein_structure_load.py @@ -16,6 +16,13 @@ from backend.protzilla.constants import paths from backend.protzilla.constants.protzilla_logging import logger +from backend.protzilla.constants.cif_columns import ( + ATOM_SITE_PREFIX, + ATOM_SITE_COLUMNS, + ATOM_SITE_COLUMNS_NUMERIC, + CHEM_COMP_PREFIX, + CHEM_COMP_COLUMNS, +) from backend.protzilla.importing.fasta_import import fasta_import from backend.protzilla.networking import download_file_from_url from backend.protzilla.utilities.utilities import copy_file_to_directory @@ -108,26 +115,54 @@ def read_alphafold_mmcif(path: Path) -> pd.DataFrame: block = doc.sole_block() - cat_name = "_atom_site." - if cat_name not in block.get_mmcif_category_names(): + if ATOM_SITE_PREFIX not in block.get_mmcif_category_names(): return pd.DataFrame() - table = block.find_mmcif_category(cat_name) - - columns = list(table.tags) - nrows = len(table) - data = {} - for j, col in enumerate(columns): - col_values = [] - for i in range(nrows): - row = table[i] - if j < len(row): - col_values.append(row[j]) - else: - col_values.append(None) - data[col] = col_values + atom_site_table = block.find_mmcif_category(ATOM_SITE_PREFIX) + + atom_site_df = pd.DataFrame( + list(atom_site_table), + columns=list(atom_site_table.tags), + dtype=pd.StringDtype(), + ) + + # convert to numeric dtype for numeric columns present in the dataframe + present_numeric_columns = [ + column for column in ATOM_SITE_COLUMNS_NUMERIC if column in atom_site_table.tags + ] + atom_site_df[present_numeric_columns] = atom_site_df[present_numeric_columns].apply( + pd.to_numeric, errors="coerce" + ) + + atom_site_df = atom_site_df.convert_dtypes() + + if CHEM_COMP_PREFIX not in block.get_mmcif_category_names(): + raise ValueError( + f"Required table with prefix {CHEM_COMP_PREFIX} not found in {path}" + ) + + chem_comp_table = block.find_mmcif_category(CHEM_COMP_PREFIX) + + chem_comp_df = pd.DataFrame( + list(chem_comp_table), + columns=list(chem_comp_table.tags), + dtype=pd.StringDtype(), + )[[CHEM_COMP_COLUMNS.ID, CHEM_COMP_COLUMNS.MON_NSTD_FLAG]] + + # convert flags to native booleans + bool_map = {"y": True, "n": False, ".": pd.NA} + + chem_comp_df[CHEM_COMP_COLUMNS.MON_NSTD_FLAG] = ( + chem_comp_df[CHEM_COMP_COLUMNS.MON_NSTD_FLAG].map(bool_map).astype("boolean") + ) - return pd.DataFrame(data) + # merge on the comp_id and drop the duplicate column + return atom_site_df.merge( + chem_comp_df, + how="left", + left_on=ATOM_SITE_COLUMNS.LABEL_COMP_ID, + right_on=CHEM_COMP_COLUMNS.ID, + ).drop(CHEM_COMP_COLUMNS.ID, axis=1) def get_correct_af_directories( diff --git a/backend/tests/protzilla/importing/test_alphafold_protein_structure_load.py b/backend/tests/protzilla/importing/test_alphafold_protein_structure_load.py index 940bbcea2..b32fe139b 100644 --- a/backend/tests/protzilla/importing/test_alphafold_protein_structure_load.py +++ b/backend/tests/protzilla/importing/test_alphafold_protein_structure_load.py @@ -29,6 +29,12 @@ check_success_of_get_df, ) from backend.protzilla.constants import paths +from backend.protzilla.constants.cif_columns import ( + ATOM_SITE_PREFIX, + ATOM_SITE_COLUMNS, + CHEM_COMP_COLUMNS, +) +from backend.protzilla.constants.data_types import DataKey def test_to_fasta_default_header_and_newline(): @@ -90,11 +96,18 @@ def test_read_alphafold_mmcif_valid_atom_site(tmp_path): """ data_test loop_ +_chem_comp.id +_chem_comp.mon_nstd_flag +SER y +# +loop_ _atom_site.id _atom_site.type_symbol +_atom_site.label_atom_id +_atom_site.label_comp_id _atom_site.Cartn_x -N N 1.0 -CA C 2.0 +1 N N SER 1.0 +2 C CA SER 2.0 """ ) @@ -102,14 +115,19 @@ def test_read_alphafold_mmcif_valid_atom_site(tmp_path): assert isinstance(df, pd.DataFrame) assert list(df.columns) == [ - "_atom_site.id", - "_atom_site.type_symbol", - "_atom_site.Cartn_x", + ATOM_SITE_COLUMNS.ID, + ATOM_SITE_COLUMNS.TYPE_SYMBOL, + ATOM_SITE_COLUMNS.LABEL_ATOM_ID, + ATOM_SITE_COLUMNS.LABEL_COMP_ID, + ATOM_SITE_COLUMNS.CARTN_X, + CHEM_COMP_COLUMNS.MON_NSTD_FLAG, ] assert len(df) == 2 - assert df["_atom_site.id"].tolist() == ["N", "CA"] - assert df["_atom_site.type_symbol"].tolist() == ["N", "C"] - assert df["_atom_site.Cartn_x"].tolist() == ["1.0", "2.0"] + assert df[ATOM_SITE_COLUMNS.ID].tolist() == [1, 2] + assert df[ATOM_SITE_COLUMNS.TYPE_SYMBOL].tolist() == ["N", "C"] + assert df[ATOM_SITE_COLUMNS.LABEL_ATOM_ID].tolist() == ["N", "CA"] + assert df[ATOM_SITE_COLUMNS.CARTN_X].tolist() == [1.0, 2.0] + assert df[CHEM_COMP_COLUMNS.MON_NSTD_FLAG].tolist() == [True, True] def test_fetch_alphafold_protein_structure_wrong_uniprot_id(): @@ -127,11 +145,11 @@ def test_fetch_alphafold_returned_keys(tmp_path, monkeypatch): out = fetch_alphafold_protein_structure("Q8WP00", persist_upload=True) assert out.keys() == { - "structure_metadata_df", - "cif_df", - "pae_df", - "plddt_df", - "amino_acid_sequences_df", + DataKey.STRUCTURE_METADATA_DF, + DataKey.CIF_DF, + DataKey.PAE_DF, + DataKey.PLDDT_DF, + DataKey.AMINO_ACID_SEQUENCES_DF, "messages", "visualization", } @@ -146,16 +164,16 @@ def test_fetch_alphafold_monomer_metadata(tmp_path, monkeypatch): ) out = fetch_alphafold_protein_structure("Q8WP00", persist_upload=True) - assert isinstance(out["structure_metadata_df"], pd.DataFrame) - assert not out["structure_metadata_df"].empty - assert out["structure_metadata_df"].iloc[0]["uniprot_accession"] == "Q8WP00" + assert isinstance(out[DataKey.STRUCTURE_METADATA_DF], pd.DataFrame) + assert not out[DataKey.STRUCTURE_METADATA_DF].empty + assert out[DataKey.STRUCTURE_METADATA_DF].iloc[0]["uniprot_accession"] == "Q8WP00" assert ( - out["structure_metadata_df"].iloc[0]["model_created_date"] + out[DataKey.STRUCTURE_METADATA_DF].iloc[0]["model_created_date"] == "2025-08-01T00:00:00Z" ) - assert out["structure_metadata_df"].iloc[0]["gene"] == "PRM1" + assert out[DataKey.STRUCTURE_METADATA_DF].iloc[0]["gene"] == "PRM1" assert ( - out["structure_metadata_df"].iloc[0]["model_used"] + out[DataKey.STRUCTURE_METADATA_DF].iloc[0]["model_used"] == "AlphaFold Monomer v2.0 pipeline" ) @@ -197,20 +215,20 @@ def test_fetch_alphafold_dfs_exist(tmp_path, monkeypatch): out = fetch_alphafold_protein_structure("Q8WP00", persist_upload=True) - cif_df = out["cif_df"] + cif_df = out[DataKey.CIF_DF] assert isinstance(cif_df, pd.DataFrame) assert not cif_df.empty - assert any(col.startswith("_atom_site.") for col in cif_df.columns) + assert any(col.startswith(ATOM_SITE_PREFIX) for col in cif_df.columns) - pae_df = out["pae_df"] + pae_df = out[DataKey.PAE_DF] assert isinstance(pae_df, pd.DataFrame) assert not pae_df.empty - plddt_df = out["plddt_df"] + plddt_df = out[DataKey.PLDDT_DF] assert isinstance(plddt_df, pd.DataFrame) assert not plddt_df.empty - seq_df = out["amino_acid_sequences_df"] + seq_df = out[DataKey.AMINO_ACID_SEQUENCES_DF] assert isinstance(seq_df, pd.DataFrame) assert not seq_df.empty @@ -280,11 +298,18 @@ def test_get_prot_structure_dfs_success(tmp_path, monkeypatch): """ data_test loop_ +_chem_comp.id +_chem_comp.mon_nstd_flag +SER y +# +loop_ _atom_site.id _atom_site.type_symbol +_atom_site.label_atom_id +_atom_site.label_comp_id _atom_site.Cartn_x -N N 1.0 -CA C 2.0 +1 N N SER 1.0 +2 C CA SER 2.0 """ ) @@ -303,34 +328,41 @@ def test_get_prot_structure_dfs_success(tmp_path, monkeypatch): out = get_monomer_structure_dfs("Q8WP00") - assert isinstance(out["structure_metadata_df"], pd.DataFrame) - assert not out["structure_metadata_df"].empty - assert out["structure_metadata_df"].iloc[0]["entry_id"] == "Q8WP00" + assert isinstance(out[DataKey.STRUCTURE_METADATA_DF], pd.DataFrame) + assert not out[DataKey.STRUCTURE_METADATA_DF].empty + assert out[DataKey.STRUCTURE_METADATA_DF].iloc[0]["entry_id"] == "Q8WP00" - assert isinstance(out["cif_df"], pd.DataFrame) - assert not out["cif_df"].empty - assert list(out["cif_df"].columns) == [ - "_atom_site.id", - "_atom_site.type_symbol", - "_atom_site.Cartn_x", + cif_df = out[DataKey.CIF_DF] + assert isinstance(cif_df, pd.DataFrame) + assert not cif_df.empty + assert list(cif_df.columns) == [ + ATOM_SITE_COLUMNS.ID, + ATOM_SITE_COLUMNS.TYPE_SYMBOL, + ATOM_SITE_COLUMNS.LABEL_ATOM_ID, + ATOM_SITE_COLUMNS.LABEL_COMP_ID, + ATOM_SITE_COLUMNS.CARTN_X, + CHEM_COMP_COLUMNS.MON_NSTD_FLAG, ] - assert out["cif_df"]["_atom_site.id"].tolist() == ["N", "CA"] - assert out["cif_df"]["_atom_site.type_symbol"].tolist() == ["N", "C"] - assert out["cif_df"]["_atom_site.Cartn_x"].tolist() == ["1.0", "2.0"] - - assert isinstance(out["pae_df"], pd.DataFrame) - assert not out["pae_df"].empty - assert out["pae_df"]["predicted_aligned_error"].tolist() == [0.1] - - assert isinstance(out["plddt_df"], pd.DataFrame) - assert not out["plddt_df"].empty - assert out["plddt_df"]["residueNumber"].tolist() == [1] - assert out["plddt_df"]["confidenceScore"].tolist() == [90] - - assert isinstance(out["amino_acid_sequences_df"], pd.DataFrame) - assert not out["amino_acid_sequences_df"].empty - assert out["amino_acid_sequences_df"]["Protein ID"].tolist() == ["Q8WP00-1"] - assert out["amino_acid_sequences_df"]["Protein Sequence"].tolist() == ["AAAA"] + assert len(cif_df) == 2 + assert cif_df[ATOM_SITE_COLUMNS.ID].tolist() == [1, 2] + assert cif_df[ATOM_SITE_COLUMNS.TYPE_SYMBOL].tolist() == ["N", "C"] + assert cif_df[ATOM_SITE_COLUMNS.LABEL_ATOM_ID].tolist() == ["N", "CA"] + assert cif_df[ATOM_SITE_COLUMNS.CARTN_X].tolist() == [1.0, 2.0] + assert cif_df[CHEM_COMP_COLUMNS.MON_NSTD_FLAG].tolist() == [True, True] + + assert isinstance(out[DataKey.PAE_DF], pd.DataFrame) + assert not out[DataKey.PAE_DF].empty + assert out[DataKey.PAE_DF]["predicted_aligned_error"].tolist() == [0.1] + + assert isinstance(out[DataKey.PLDDT_DF], pd.DataFrame) + assert not out[DataKey.PLDDT_DF].empty + assert out[DataKey.PLDDT_DF]["residueNumber"].tolist() == [1] + assert out[DataKey.PLDDT_DF]["confidenceScore"].tolist() == [90] + + assert isinstance(out[DataKey.AMINO_ACID_SEQUENCES_DF], pd.DataFrame) + assert not out[DataKey.AMINO_ACID_SEQUENCES_DF].empty + assert out[DataKey.AMINO_ACID_SEQUENCES_DF]["Protein ID"].tolist() == ["Q8WP00-1"] + assert out[DataKey.AMINO_ACID_SEQUENCES_DF]["Protein Sequence"].tolist() == ["AAAA"] assert any(d.get("level") == logging.INFO for d in out["messages"]) or any( "Successfully loaded" in d.get("msg", "") for d in out["messages"] @@ -421,11 +453,13 @@ def test_get_amino_acid_sequences_df_and_handle_files(tmp_path, monkeypatch): out = handle_alphafold_files( {}, "P", "TESTSEQ", metadata_df, "P", persist_upload=False ) - assert "amino_acid_sequences_df" in out - assert isinstance(out["cif_df"], pd.DataFrame) and out["cif_df"].empty - assert isinstance(out["pae_df"], pd.DataFrame) and out["pae_df"].empty - assert isinstance(out["plddt_df"], pd.DataFrame) and out["plddt_df"].empty - assert isinstance(out["amino_acid_sequences_df"], pd.DataFrame) + assert DataKey.AMINO_ACID_SEQUENCES_DF in out + assert isinstance(out[DataKey.CIF_DF], pd.DataFrame) and out[DataKey.CIF_DF].empty + assert isinstance(out[DataKey.PAE_DF], pd.DataFrame) and out[DataKey.PAE_DF].empty + assert ( + isinstance(out[DataKey.PLDDT_DF], pd.DataFrame) and out[DataKey.PLDDT_DF].empty + ) + assert isinstance(out[DataKey.AMINO_ACID_SEQUENCES_DF], pd.DataFrame) def test_upload_multimer_prediction_basic(tmp_path, monkeypatch): @@ -439,9 +473,16 @@ def test_upload_multimer_prediction_basic(tmp_path, monkeypatch): """ data_test loop_ +_chem_comp.id +_chem_comp.mon_nstd_flag +SER y +# +loop_ _atom_site.id _atom_site.type_symbol -N N +_atom_site.label_atom_id +_atom_site.label_comp_id +1 N N SER """ ) conf = tmp_path / "conf.json" @@ -494,39 +535,48 @@ def _copy(src, dest_dir): persist_upload=True, ) - assert isinstance(out["structure_metadata_df"], pd.DataFrame) + assert isinstance(out[DataKey.STRUCTURE_METADATA_DF], pd.DataFrame) # check metadata contents - mdf = out["structure_metadata_df"] + mdf = out[DataKey.STRUCTURE_METADATA_DF] assert mdf.iloc[0]["entry_id"] == "M1" assert mdf.iloc[0]["uniprot_ids"] == ["X"] assert mdf.iloc[0]["model_used"] == "m" # cif contents - cif_df = out["cif_df"] + cif_df = out[DataKey.CIF_DF] assert isinstance(cif_df, pd.DataFrame) - assert list(cif_df.columns) == ["_atom_site.id", "_atom_site.type_symbol"] - assert cif_df["_atom_site.id"].tolist() == ["N"] - assert cif_df["_atom_site.type_symbol"].tolist() == ["N"] + assert list(cif_df.columns) == [ + ATOM_SITE_COLUMNS.ID, + ATOM_SITE_COLUMNS.TYPE_SYMBOL, + ATOM_SITE_COLUMNS.LABEL_ATOM_ID, + ATOM_SITE_COLUMNS.LABEL_COMP_ID, + CHEM_COMP_COLUMNS.MON_NSTD_FLAG, + ] + assert cif_df[ATOM_SITE_COLUMNS.ID].tolist() == [1] + assert cif_df[ATOM_SITE_COLUMNS.TYPE_SYMBOL].tolist() == ["N"] + assert cif_df[ATOM_SITE_COLUMNS.LABEL_ATOM_ID].tolist() == ["N"] + assert cif_df[ATOM_SITE_COLUMNS.LABEL_COMP_ID].tolist() == ["SER"] + assert cif_df[CHEM_COMP_COLUMNS.MON_NSTD_FLAG].tolist() == [True] # confidence JSON - conf_df = out["confidence_df"] + conf_df = out[DataKey.CONFIDENCE_DF] assert isinstance(conf_df, pd.DataFrame) assert conf_df["residueNumber"].tolist() == [1] assert conf_df["confidenceScore"].tolist() == [99] # full data normalization - full_df = out["full_data_df"] + full_df = out[DataKey.FULL_DATA_DF] assert isinstance(full_df, pd.DataFrame) assert full_df.iloc[0]["a"] == [1, 2] # job request JSON - job_df = out["job_request_df"] + job_df = out[DataKey.JOB_REQUEST_DF] assert isinstance(job_df, pd.DataFrame) assert job_df.iloc[0]["name"] == "test_job" assert job_df.iloc[0]["dialect"] == "alphafoldserver" # sequences - seqs = out["amino_acid_sequences_df"] + seqs = out[DataKey.AMINO_ACID_SEQUENCES_DF] assert isinstance(seqs, pd.DataFrame) assert seqs["Protein Sequence"].tolist() == ["AAAA"] assert any(str(v).startswith("X") for v in seqs["Protein ID"].tolist()) @@ -612,7 +662,9 @@ def test_upload_multimer_prediction_no_persist(tmp_path, monkeypatch): fasta = tmp_path / "seqs.fasta" fasta.write_text(">alpha|X\nAAAA\n") cif = tmp_path / "m.cif" - cif.write_text("data_test\nloop_\n_atom_site.id\nN\n") + cif.write_text( + "data_test\nloop_\n_chem_comp.id\n_chem_comp.mon_nstd_flag\nSER y\nloop_\n#\n_atom_site.id\n_atom_site.label_comp_id\nN SER\n" + ) conf = tmp_path / "conf.json" conf.write_text('[{"residueNumber":1, "confidenceScore":99}]') full = tmp_path / "full.json" @@ -653,10 +705,10 @@ def test_upload_multimer_prediction_no_persist(tmp_path, monkeypatch): ) # verify dataframes are returned - assert isinstance(out["structure_metadata_df"], pd.DataFrame) - assert isinstance(out["cif_df"], pd.DataFrame) - assert isinstance(out["job_request_df"], pd.DataFrame) - assert out["job_request_df"].iloc[0]["name"] == "test_job_2" + assert isinstance(out[DataKey.STRUCTURE_METADATA_DF], pd.DataFrame) + assert isinstance(out[DataKey.CIF_DF], pd.DataFrame) + assert isinstance(out[DataKey.JOB_REQUEST_DF], pd.DataFrame) + assert out[DataKey.JOB_REQUEST_DF].iloc[0]["name"] == "test_job_2" # directory should still exist (created for the entry) upload_dir = tmp_path / "M2" assert not upload_dir.exists() @@ -692,7 +744,9 @@ def test_get_prot_structure_dfs_missing_fasta(tmp_path, monkeypatch): # create CIF but no FASTA cif = prot_dir / "test.cif" - cif.write_text("data_test\nloop_\n_atom_site.id\nN\n") + cif.write_text( + "data_test\nloop_\n_chem_comp.id\n_chem_comp.mon_nstd_flag\nSER y\nloop_\n#\n_atom_site.id\n_atom_site.label_comp_id\nN SER\n" + ) with pytest.raises(FileNotFoundError, match="No FASTA file found"): get_monomer_structure_dfs("NOFASTA") @@ -712,7 +766,9 @@ def test_get_prot_structure_dfs_missing_json(tmp_path, monkeypatch): # create CIF and FASTA but no JSON cif = prot_dir / "test.cif" - cif.write_text("data_test\nloop_\n_atom_site.id\nN\n") + cif.write_text( + "data_test\nloop_\n_chem_comp.id\n_chem_comp.mon_nstd_flag\nSER y\nloop_\n#\n_atom_site.id\n_atom_site.label_comp_id\nN SER\n" + ) fasta = prot_dir / "test.fasta" # valid header for parse_fasta_id (expects at least one "|" in the id) @@ -820,18 +876,30 @@ def test_get_cif_df_from_disk_multiple_cif_warns(tmp_path): """ data_test loop_ +_chem_comp.id +_chem_comp.mon_nstd_flag +SER y +# +loop_ _atom_site.id _atom_site.type_symbol -N N +_atom_site.label_comp_id +N N SER """ ) cif2.write_text( """ data_test loop_ +_chem_comp.id +_chem_comp.mon_nstd_flag +SER y +# +loop_ _atom_site.id _atom_site.type_symbol -CA C +_atom_site.label_comp_id +CA C SER """ ) @@ -872,9 +940,15 @@ def test_get_multimer_structure_dfs_success(tmp_path, monkeypatch): """ data_test loop_ +_chem_comp.id +_chem_comp.mon_nstd_flag +SER y +# +loop_ _atom_site.id _atom_site.type_symbol -N N +_atom_site.label_comp_id +N N SER """ ) @@ -909,17 +983,17 @@ def test_get_multimer_structure_dfs_success(tmp_path, monkeypatch): ) out = get_multimer_structure_dfs("M1") - assert isinstance(out["structure_metadata_df"], pd.DataFrame) - assert isinstance(out["cif_df"], pd.DataFrame) - assert isinstance(out["amino_acid_sequences_df"], pd.DataFrame) - assert isinstance(out["confidence_df"], pd.DataFrame) - assert isinstance(out["full_data_df"], pd.DataFrame) - assert isinstance(out["job_request_df"], pd.DataFrame) - - assert "chain_iptm" in out["confidence_df"].columns - assert "pae" in out["full_data_df"].columns - assert out["job_request_df"].iloc[0]["name"] == "multimer_job" - assert out["job_request_df"].iloc[0]["version"] == 3 + assert isinstance(out[DataKey.STRUCTURE_METADATA_DF], pd.DataFrame) + assert isinstance(out[DataKey.CIF_DF], pd.DataFrame) + assert isinstance(out[DataKey.AMINO_ACID_SEQUENCES_DF], pd.DataFrame) + assert isinstance(out[DataKey.CONFIDENCE_DF], pd.DataFrame) + assert isinstance(out[DataKey.FULL_DATA_DF], pd.DataFrame) + assert isinstance(out[DataKey.JOB_REQUEST_DF], pd.DataFrame) + + assert "chain_iptm" in out[DataKey.CONFIDENCE_DF].columns + assert "pae" in out[DataKey.FULL_DATA_DF].columns + assert out[DataKey.JOB_REQUEST_DF].iloc[0]["name"] == "multimer_job" + assert out[DataKey.JOB_REQUEST_DF].iloc[0]["version"] == 3 assert any(m.get("level") == logging.INFO for m in out["messages"]) or any( "Successfully loaded" in str(m.get("msg", "")) for m in out["messages"] @@ -956,9 +1030,15 @@ def test_get_multimer_structure_dfs_json_fallback_warns(tmp_path, monkeypatch): """ data_test loop_ +_chem_comp.id +_chem_comp.mon_nstd_flag +SER y +# +loop_ _atom_site.id _atom_site.type_symbol -N N +_atom_site.label_comp_id +N N SER """ ) From 27a3c867dd27414977786317f7225645cac1301c Mon Sep 17 00:00:00 2001 From: 3dot141592 Date: Mon, 11 May 2026 13:39:53 +0200 Subject: [PATCH 06/22] Add descriptions and rename --- backend/protzilla/data_analysis/amino_acid_spheres.py | 4 ++-- backend/protzilla/data_analysis/geometry_operations.py | 8 +++++++- 2 files changed, 9 insertions(+), 3 deletions(-) diff --git a/backend/protzilla/data_analysis/amino_acid_spheres.py b/backend/protzilla/data_analysis/amino_acid_spheres.py index c9cd57ace..ad7d00527 100644 --- a/backend/protzilla/data_analysis/amino_acid_spheres.py +++ b/backend/protzilla/data_analysis/amino_acid_spheres.py @@ -5,7 +5,7 @@ import trimesh from backend.protzilla.data_analysis.geometry_operations import ( - calculate_centroid, + calculate_center_point, extract_points_from_cif, find_farthest_point, mesh_to_polyhedron, @@ -34,7 +34,7 @@ def calculate_amino_acid_spheres( residue_range=(residue_position, residue_position), chain_id=chain_id, ) - center = calculate_centroid(residue_points) + center = calculate_center_point(residue_points) furthest_point = find_farthest_point(residue_points, center) radius = float(np.linalg.norm(furthest_point - center)) diff --git a/backend/protzilla/data_analysis/geometry_operations.py b/backend/protzilla/data_analysis/geometry_operations.py index 7b22e29cd..b85d35315 100644 --- a/backend/protzilla/data_analysis/geometry_operations.py +++ b/backend/protzilla/data_analysis/geometry_operations.py @@ -118,7 +118,10 @@ def point_cloud_to_polyhedron(points: np.ndarray) -> dict: return mesh_to_polyhedron(build_convex_hull(points)) -def calculate_centroid(points: np.ndarray) -> np.ndarray: +def calculate_center_point(points: np.ndarray) -> np.ndarray: + """ + Calculate the center point of a point cloud. + """ if points.ndim != 2 or points.shape[1] != 3: raise ValueError(f"Expected points with shape (n, 3), got {points.shape}.") if len(points) == 0: @@ -128,6 +131,9 @@ def calculate_centroid(points: np.ndarray) -> np.ndarray: def find_farthest_point(points: np.ndarray, reference_point: np.ndarray) -> np.ndarray: + """ + Calculate the farthest point inside a point cloud from a reference point. + """ if points.ndim != 2 or points.shape[1] != 3: raise ValueError(f"Expected points with shape (n, 3), got {points.shape}.") if len(points) == 0: From 84857a5cc47a63480903282857921429ceb66e1e Mon Sep 17 00:00:00 2001 From: 3dot141592 Date: Mon, 11 May 2026 14:21:03 +0200 Subject: [PATCH 07/22] WIP: fdw distance --- backend/protzilla/constants/option_types.py | 1 + backend/protzilla/constants/van_der_waals.py | 124 ++++++++++++++++++ .../data_analysis/geometry_operations.py | 16 +++ 3 files changed, 141 insertions(+) create mode 100644 backend/protzilla/constants/van_der_waals.py diff --git a/backend/protzilla/constants/option_types.py b/backend/protzilla/constants/option_types.py index 678efe9cb..8a446fae3 100644 --- a/backend/protzilla/constants/option_types.py +++ b/backend/protzilla/constants/option_types.py @@ -14,6 +14,7 @@ class LogBaseWithNoneType(StrEnum): LOG2 = "log2" LOG10 = "log10" # i hate this ~T + # lmao ~R NONE = "None" diff --git a/backend/protzilla/constants/van_der_waals.py b/backend/protzilla/constants/van_der_waals.py new file mode 100644 index 000000000..a7e25ad8c --- /dev/null +++ b/backend/protzilla/constants/van_der_waals.py @@ -0,0 +1,124 @@ +vdw_radii = { + """ + Van der Waals radius in Ångström + """ + + "H": 1.20, + "He": 1.43, + "Li": 2.12, + "Be": 1.98, + "B": 1.91, + "C": 1.77, + "N": 1.66, + "O": 1.50, + "F": 1.46, + "Ne": 1.58, + "Na": 2.50, + "Mg": 2.51, + "Al": 2.25, + "Si": 2.19, + "P": 1.90, + "S": 1.89, + "Cl": 1.82, + "Ar": 1.83, + "K": 2.73, + "Ca": 2.62, + "Sc": 2.58, + "Ti": 2.46, + "V": 2.42, + "Cr": 2.45, + "Mn": 2.45, + "Fe": 2.44, + "Co": 2.40, + "Ni": 2.40, + "Cu": 2.38, + "Zn": 2.39, + "Ga": 2.32, + "Ge": 2.29, + "As": 1.88, + "Se": 1.82, + "Br": 1.86, + "Kr": 2.25, + "Rb": 3.21, + "Sr": 2.84, + "Y": 2.75, + "Zr": 2.52, + "Nb": 2.56, + "Mo": 2.45, + "Tc": 2.44, + "Ru": 2.46, + "Rh": 2.44, + "Pd": 2.15, + "Ag": 2.53, + "Cd": 2.49, + "In": 2.43, + "Sn": 2.42, + "Sb": 2.47, + "Te": 1.99, + "I": 2.04, + "Xe": 2.06, + "Cs": 3.48, + "Ba": 3.03, + "La": 2.98, + "Ce": 2.88, + "Pr": 2.92, + "Nd": 2.95, + "Pm": None, + "Sm": 2.90, + "Eu": 2.87, + "Gd": 2.83, + "Tb": 2.79, + "Dy": 2.87, + "Ho": 2.81, + "Er": 2.83, + "Tm": 2.79, + "Yb": 2.80, + "Lu": 2.74, + "Hf": 2.63, + "Ta": 2.53, + "W": 2.57, + "Re": 2.49, + "Os": 2.48, + "Ir": 2.41, + "Pt": 2.29, + "Au": 2.32, + "Hg": 2.45, + "Tl": 2.47, + "Pb": 2.60, + "Bi": 2.54, + "Po": None, + "At": None, + "Rn": None, + "Fr": None, + "Ra": None, + "Ac": 2.80, + "Th": 2.93, + "Pa": 2.88, + "U": 2.71, + "Np": 2.82, + "Pu": 2.81, + "Am": 2.83, + "Cm": 3.05, + "Bk": 3.40, + "Cf": 3.05, + "Es": 2.70, + "Fm": None, + "Md": None, + "No": None, + "Lr": None, + "Rf": None, + "Db": None, + "Sg": None, + "Bh": None, + "Hs": None, + "Mt": None, + "Ds": None, + "Rg": None, + "Cn": None, + "Nh": None, + "Fl": None, + "Mc": None, + "Lv": None, + "Ts": None, + "Og": None, +} \ No newline at end of file diff --git a/backend/protzilla/data_analysis/geometry_operations.py b/backend/protzilla/data_analysis/geometry_operations.py index b85d35315..9605985e2 100644 --- a/backend/protzilla/data_analysis/geometry_operations.py +++ b/backend/protzilla/data_analysis/geometry_operations.py @@ -146,6 +146,22 @@ def find_farthest_point(points: np.ndarray, reference_point: np.ndarray) -> np.n distances = np.linalg.norm(points - reference_point, axis=1) return points[np.argmax(distances)] +def find_farthest_point_fdw(points: np.ndarray, reference_point: np.ndarray) -> np.ndarray: + """ + Calculate the point, whose fdw-"Bubble" is most distant to a reference point. + """ + if points.ndim != 2 or points.shape[1] != 3: + raise ValueError(f"Expected points with shape (n, 3), got {points.shape}.") + if len(points) == 0: + raise ValueError("At least one point is required to calculate a maximum distance.") + if reference_point.shape != (3,): + raise ValueError( + f"Expected reference_point with shape (3,), got {reference_point.shape}." + ) + + distances = np.linalg.norm(points - reference_point, axis=1) + return points[np.argmax(distances)] + def meshes_intersect( mesh_a: trimesh.Trimesh, From 9cf5e4128d05f6200a068a8e6982cffd801f95cf Mon Sep 17 00:00:00 2001 From: Elena-kal Date: Mon, 11 May 2026 14:40:26 +0200 Subject: [PATCH 08/22] Add default step for ptm validation --- backend/protzilla/all_steps.py | 1 + backend/protzilla/methods/data_analysis.py | 22 ++++++++++++++++------ 2 files changed, 17 insertions(+), 6 deletions(-) diff --git a/backend/protzilla/all_steps.py b/backend/protzilla/all_steps.py index 7392a6bd4..98a868492 100644 --- a/backend/protzilla/all_steps.py +++ b/backend/protzilla/all_steps.py @@ -82,6 +82,7 @@ data_analysis.PTMDetailsVisualization, data_analysis.CrosslinkingValidationWithAngstromDeviation, data_analysis.CrosslinkingValidationWithAngstromDeviationForMultimer, + data_analysis.PtmValidation, data_preprocessing.ImputationByMinPerSample, data_integration.EnrichmentAnalysisGOAnalysisWithString, data_integration.EnrichmentAnalysisGOAnalysisWithEnrichr, diff --git a/backend/protzilla/methods/data_analysis.py b/backend/protzilla/methods/data_analysis.py index 51d487642..1872c6d34 100644 --- a/backend/protzilla/methods/data_analysis.py +++ b/backend/protzilla/methods/data_analysis.py @@ -95,12 +95,8 @@ monomer_validation, multimer_validation, ) -from backend.protzilla.run import Run -from backend.protzilla.methods.importing import ( - ImportMonomerStructurePredictionFromDisk, - AlphaFoldPredictionLoad, - ImportMultimerStructurePredictionFromDisk, - UploadMultimerPredictions, +from protzilla.data_analysis.ptm_analysis import( + ptm_validation ) @@ -2452,3 +2448,17 @@ def create_form(self): label="Ångström Deviation - Multimer", input_fields=[], ) + + +class PtmValidation(DataAnalysisStep): + display_name = "PTM Validation" + operation = "Peptide analysis" + method_description = "Validates PTMs in protein structure predictions." + output_keys = ["ptm_validation_df"] + calc_method = staticmethod(ptm_validation) + + def create_form(self): + return Form( + label="PTM Validation", + input_fields=[], + ) From 7be7535dc2fae9b1ac8f8871e5ba6ef6c1e9f5ea Mon Sep 17 00:00:00 2001 From: Elena-kal Date: Mon, 11 May 2026 14:41:25 +0200 Subject: [PATCH 09/22] add functions to obtain ptms with all their atoms and their positions and their centerpoints and radius --- .../protzilla/data_analysis/ptm_analysis.py | 7 ++ backend/protzilla/utilities/ptm_helper.py | 96 +++++++++++++++++++ 2 files changed, 103 insertions(+) create mode 100644 backend/protzilla/utilities/ptm_helper.py diff --git a/backend/protzilla/data_analysis/ptm_analysis.py b/backend/protzilla/data_analysis/ptm_analysis.py index 7174e8a09..2c170eab3 100644 --- a/backend/protzilla/data_analysis/ptm_analysis.py +++ b/backend/protzilla/data_analysis/ptm_analysis.py @@ -3,6 +3,7 @@ import re from backend.protzilla.utilities.transform_dfs import long_to_wide +from backend.protzilla.utilities.ptm_helper import get_all_ptm_atoms_with_coordinates, get_center_points_and_radius_for_each_ptm def ptms_per_sample(psm_df: pd.DataFrame) -> dict: @@ -91,3 +92,9 @@ def from_string(mod_string: str) -> tuple[int, str]: name = name[1:] if name[0] == " " else name return amount, name + + +def ptm_validation(cif_df: pd.DataFrame): + ptm_list = get_all_ptm_atoms_with_coordinates(cif_df) + ptm_list = get_center_points_and_radius_for_each_ptm(ptm_list) + return {"ptm_validation_df": pd.DataFrame(ptm_list)} \ No newline at end of file diff --git a/backend/protzilla/utilities/ptm_helper.py b/backend/protzilla/utilities/ptm_helper.py new file mode 100644 index 000000000..3ad41ae20 --- /dev/null +++ b/backend/protzilla/utilities/ptm_helper.py @@ -0,0 +1,96 @@ +import pandas as pd +import numpy as np + + +def get_all_ptm_atoms_with_coordinates(cif_df: pd.DataFrame) -> list: + """ + Extracts all non-standard residues (PTMs) and their atomic coordinates + from a parsed mmCIF DataFrame. + """ + ptm_list = [] + + is_non_standard = cif_df["_chem_comp.mon_nstd_flag"] == False + ptm_atoms_df = cif_df[is_non_standard] + + grouped_ptms = ptm_atoms_df.groupby( + ["_atom_site.label_asym_id", "_atom_site.label_seq_id"] + ) + + for (chain_id, seq_id), atom_group in grouped_ptms: + comp_id = atom_group["_atom_site.label_comp_id"].iloc[0] + atoms_subset = atom_group[ + [ + "_atom_site.label_atom_id", + "_atom_site.Cartn_x", + "_atom_site.Cartn_y", + "_atom_site.Cartn_z", + ] + ].rename( + columns={ + "_atom_site.label_atom_id": "atom_name", + "_atom_site.Cartn_x": "x", + "_atom_site.Cartn_y": "y", + "_atom_site.Cartn_z": "z", + } + ) + + atoms_data = atoms_subset.to_dict(orient="records") + + ptm_list.append( + { + "ptm_name": comp_id, + "chain": chain_id, + "position": seq_id, + "atoms": atoms_data, + } + ) + + return ptm_list + + +def calculate_center_point(points: np.ndarray) -> np.ndarray: + """ + Calculate the center point of a point cloud. + """ + if points.ndim != 2 or points.shape[1] != 3: + raise ValueError(f"Expected points with shape (n, 3), got {points.shape}.") + if len(points) == 0: + raise ValueError("At least one point is required to calculate a centroid.") + + return points.mean(axis=0) + + +def find_farthest_point(points: np.ndarray, reference_point: np.ndarray) -> np.ndarray: + """ + Calculate the farthest point inside a point cloud from a reference point. + """ + if points.ndim != 2 or points.shape[1] != 3: + raise ValueError(f"Expected points with shape (n, 3), got {points.shape}.") + if len(points) == 0: + raise ValueError("At least one point is required to calculate a maximum distance.") + if reference_point.shape != (3,): + raise ValueError( + f"Expected reference_point with shape (3,), got {reference_point.shape}." + ) + + distances = np.linalg.norm(points - reference_point, axis=1) + return points[np.argmax(distances)] + +def get_center_points_and_radius_for_each_ptm(ptm_list: list) -> list: + for ptm in ptm_list: + coords = [ + [atom["x"], atom["y"], atom["z"]] + for atom in ptm["atoms"] + ] + + coords_array = np.array(coords, dtype=np.float32) + + center_point = calculate_center_point(coords_array) + farthest_point = find_farthest_point(coords_array, center_point) + + radius = np.linalg.norm(farthest_point - center_point) + + ptm["center_point"] = center_point.tolist() + ptm["radius"] = float(radius) + + return ptm_list \ No newline at end of file From ec7e47693fb649cb4223fe0a99c24df9cce66ec0 Mon Sep 17 00:00:00 2001 From: Elena-kal Date: Mon, 11 May 2026 14:43:59 +0200 Subject: [PATCH 10/22] add tests for the ptm helper functions --- .../protzilla/utilities/test_ptm_helper.py | 98 +++++++++++++++++++ 1 file changed, 98 insertions(+) create mode 100644 backend/tests/protzilla/utilities/test_ptm_helper.py diff --git a/backend/tests/protzilla/utilities/test_ptm_helper.py b/backend/tests/protzilla/utilities/test_ptm_helper.py new file mode 100644 index 000000000..321b2e010 --- /dev/null +++ b/backend/tests/protzilla/utilities/test_ptm_helper.py @@ -0,0 +1,98 @@ +import pandas as pd +from backend.protzilla.utilities.ptm_helper import get_all_ptm_atoms_with_coordinates, get_center_points_and_radius_for_each_ptm +import pytest +import numpy as np + +def test_get_all_ptm_atoms_with_coordinates(): + data = { + "_chem_comp.mon_nstd_flag": [True, True, False, False], + "_atom_site.label_asym_id": ["A", "A", "B", "B"], + "_atom_site.label_seq_id": [10, 10, 42, 42], + "_atom_site.label_comp_id": ["MET", "MET", "SEP", "SEP"], + "_atom_site.label_atom_id": ["N", "CA", "P", "O1P"], + "_atom_site.Cartn_x": [1.0, 2.0, 10.0, 11.0], + "_atom_site.Cartn_y": [1.1, 2.1, 10.1, 11.1], + "_atom_site.Cartn_z": [1.2, 2.2, 10.2, 11.2], + } + cif_df = pd.DataFrame(data) + + result = get_all_ptm_atoms_with_coordinates(cif_df) + + assert len(result) == 1 + + ptm = result[0] + assert ptm["ptm_name"] == "SEP" + assert ptm["chain"] == "B" + assert ptm["position"] == 42 + + atoms = ptm["atoms"] + assert len(atoms) == 2 + + assert atoms[0] == { + "atom_name": "P", + "x": 10.0, + "y": 10.1, + "z": 10.2 + } + + assert atoms[1] == { + "atom_name": "O1P", + "x": 11.0, + "y": 11.1, + "z": 11.2 + } + + +def test_get_center_points_and_radius_for_each_ptm(): + ptm_list = [ + { + "ptm_name": "SYM", + "chain": "A", + "position": 1, + "atoms": [ + {"atom_name": "A1", "x": 1.0, "y": 0.0, "z": 0.0}, + {"atom_name": "A2", "x": -1.0, "y": 0.0, "z": 0.0}, + {"atom_name": "A3", "x": 0.0, "y": 1.0, "z": 0.0}, + {"atom_name": "A4", "x": 0.0, "y": -1.0, "z": 0.0}, + ] + }, + { + "ptm_name": "SGL", + "chain": "B", + "position": 2, + "atoms": [ + {"atom_name": "A1", "x": 5.0, "y": 5.0, "z": 5.0}, + ] + }, + { + "ptm_name": "AB", + "chain": "C", + "position": 3, + "atoms": [ + {"atom_name": "A1", "x": 12.0, "y": 10.0, "z": 8.0}, + {"atom_name": "A2", "x": 8.0, "y": 10.0, "z": 7.0}, + {"atom_name": "A3", "x": 10.0, "y": 10.0, "z": 15.0}, + ] + } + ] + + result_list = get_center_points_and_radius_for_each_ptm(ptm_list) + + ptm_1 = result_list[0] + + result_list = get_center_points_and_radius_for_each_ptm(ptm_list) + + ptm_1 = result_list[0] + np.testing.assert_allclose(ptm_1["center_point"], [0.0, 0.0, 0.0], atol=1e-6) + assert pytest.approx(ptm_1["radius"], 1e-6) == 1.0 + + ptm_2 = result_list[1] + np.testing.assert_allclose(ptm_2["center_point"], [5.0, 5.0, 5.0], atol=1e-6) + assert pytest.approx(ptm_2["radius"], 1e-6) == 0.0 + + ptm_3 = result_list[2] + + expected_center = [10.0, 10.0, 10.0] + np.testing.assert_allclose(ptm_3["center_point"], expected_center, atol=1e-6) + + assert pytest.approx(ptm_3["radius"], 1e-6) == 5.0 \ No newline at end of file From c211ff9459f5307f429230a4f2ca6dc767994c6d Mon Sep 17 00:00:00 2001 From: Elena-kal Date: Mon, 11 May 2026 14:44:43 +0200 Subject: [PATCH 11/22] format with black --- .../protzilla/data_analysis/ptm_analysis.py | 7 ++- backend/protzilla/methods/data_analysis.py | 4 +- backend/protzilla/utilities/ptm_helper.py | 30 +++++------ .../protzilla/utilities/test_ptm_helper.py | 50 ++++++++----------- 4 files changed, 43 insertions(+), 48 deletions(-) diff --git a/backend/protzilla/data_analysis/ptm_analysis.py b/backend/protzilla/data_analysis/ptm_analysis.py index 2c170eab3..8665eaaea 100644 --- a/backend/protzilla/data_analysis/ptm_analysis.py +++ b/backend/protzilla/data_analysis/ptm_analysis.py @@ -3,7 +3,10 @@ import re from backend.protzilla.utilities.transform_dfs import long_to_wide -from backend.protzilla.utilities.ptm_helper import get_all_ptm_atoms_with_coordinates, get_center_points_and_radius_for_each_ptm +from backend.protzilla.utilities.ptm_helper import ( + get_all_ptm_atoms_with_coordinates, + get_center_points_and_radius_for_each_ptm, +) def ptms_per_sample(psm_df: pd.DataFrame) -> dict: @@ -97,4 +100,4 @@ def from_string(mod_string: str) -> tuple[int, str]: def ptm_validation(cif_df: pd.DataFrame): ptm_list = get_all_ptm_atoms_with_coordinates(cif_df) ptm_list = get_center_points_and_radius_for_each_ptm(ptm_list) - return {"ptm_validation_df": pd.DataFrame(ptm_list)} \ No newline at end of file + return {"ptm_validation_df": pd.DataFrame(ptm_list)} diff --git a/backend/protzilla/methods/data_analysis.py b/backend/protzilla/methods/data_analysis.py index 1872c6d34..e8f334b13 100644 --- a/backend/protzilla/methods/data_analysis.py +++ b/backend/protzilla/methods/data_analysis.py @@ -95,9 +95,7 @@ monomer_validation, multimer_validation, ) -from protzilla.data_analysis.ptm_analysis import( - ptm_validation -) +from protzilla.data_analysis.ptm_analysis import ptm_validation class TTestType(Enum): diff --git a/backend/protzilla/utilities/ptm_helper.py b/backend/protzilla/utilities/ptm_helper.py index 3ad41ae20..65d6a1229 100644 --- a/backend/protzilla/utilities/ptm_helper.py +++ b/backend/protzilla/utilities/ptm_helper.py @@ -50,7 +50,7 @@ def get_all_ptm_atoms_with_coordinates(cif_df: pd.DataFrame) -> list: def calculate_center_point(points: np.ndarray) -> np.ndarray: """ - Calculate the center point of a point cloud. + Calculate the center point of a point cloud. """ if points.ndim != 2 or points.shape[1] != 3: raise ValueError(f"Expected points with shape (n, 3), got {points.shape}.") @@ -62,12 +62,14 @@ def calculate_center_point(points: np.ndarray) -> np.ndarray: def find_farthest_point(points: np.ndarray, reference_point: np.ndarray) -> np.ndarray: """ - Calculate the farthest point inside a point cloud from a reference point. + Calculate the farthest point inside a point cloud from a reference point. """ if points.ndim != 2 or points.shape[1] != 3: raise ValueError(f"Expected points with shape (n, 3), got {points.shape}.") if len(points) == 0: - raise ValueError("At least one point is required to calculate a maximum distance.") + raise ValueError( + "At least one point is required to calculate a maximum distance." + ) if reference_point.shape != (3,): raise ValueError( f"Expected reference_point with shape (3,), got {reference_point.shape}." @@ -76,21 +78,19 @@ def find_farthest_point(points: np.ndarray, reference_point: np.ndarray) -> np.n distances = np.linalg.norm(points - reference_point, axis=1) return points[np.argmax(distances)] + def get_center_points_and_radius_for_each_ptm(ptm_list: list) -> list: for ptm in ptm_list: - coords = [ - [atom["x"], atom["y"], atom["z"]] - for atom in ptm["atoms"] - ] - + coords = [[atom["x"], atom["y"], atom["z"]] for atom in ptm["atoms"]] + coords_array = np.array(coords, dtype=np.float32) - + center_point = calculate_center_point(coords_array) farthest_point = find_farthest_point(coords_array, center_point) - + radius = np.linalg.norm(farthest_point - center_point) - - ptm["center_point"] = center_point.tolist() - ptm["radius"] = float(radius) - - return ptm_list \ No newline at end of file + + ptm["center_point"] = center_point.tolist() + ptm["radius"] = float(radius) + + return ptm_list diff --git a/backend/tests/protzilla/utilities/test_ptm_helper.py b/backend/tests/protzilla/utilities/test_ptm_helper.py index 321b2e010..f63e97461 100644 --- a/backend/tests/protzilla/utilities/test_ptm_helper.py +++ b/backend/tests/protzilla/utilities/test_ptm_helper.py @@ -1,8 +1,12 @@ import pandas as pd -from backend.protzilla.utilities.ptm_helper import get_all_ptm_atoms_with_coordinates, get_center_points_and_radius_for_each_ptm +from backend.protzilla.utilities.ptm_helper import ( + get_all_ptm_atoms_with_coordinates, + get_center_points_and_radius_for_each_ptm, +) import pytest import numpy as np + def test_get_all_ptm_atoms_with_coordinates(): data = { "_chem_comp.mon_nstd_flag": [True, True, False, False], @@ -19,34 +23,24 @@ def test_get_all_ptm_atoms_with_coordinates(): result = get_all_ptm_atoms_with_coordinates(cif_df) assert len(result) == 1 - + ptm = result[0] assert ptm["ptm_name"] == "SEP" assert ptm["chain"] == "B" assert ptm["position"] == 42 - + atoms = ptm["atoms"] assert len(atoms) == 2 - - assert atoms[0] == { - "atom_name": "P", - "x": 10.0, - "y": 10.1, - "z": 10.2 - } - - assert atoms[1] == { - "atom_name": "O1P", - "x": 11.0, - "y": 11.1, - "z": 11.2 - } + + assert atoms[0] == {"atom_name": "P", "x": 10.0, "y": 10.1, "z": 10.2} + + assert atoms[1] == {"atom_name": "O1P", "x": 11.0, "y": 11.1, "z": 11.2} def test_get_center_points_and_radius_for_each_ptm(): ptm_list = [ { - "ptm_name": "SYM", + "ptm_name": "SYM", "chain": "A", "position": 1, "atoms": [ @@ -54,32 +48,32 @@ def test_get_center_points_and_radius_for_each_ptm(): {"atom_name": "A2", "x": -1.0, "y": 0.0, "z": 0.0}, {"atom_name": "A3", "x": 0.0, "y": 1.0, "z": 0.0}, {"atom_name": "A4", "x": 0.0, "y": -1.0, "z": 0.0}, - ] + ], }, { - "ptm_name": "SGL", + "ptm_name": "SGL", "chain": "B", "position": 2, "atoms": [ {"atom_name": "A1", "x": 5.0, "y": 5.0, "z": 5.0}, - ] + ], }, { - "ptm_name": "AB", + "ptm_name": "AB", "chain": "C", "position": 3, "atoms": [ {"atom_name": "A1", "x": 12.0, "y": 10.0, "z": 8.0}, {"atom_name": "A2", "x": 8.0, "y": 10.0, "z": 7.0}, {"atom_name": "A3", "x": 10.0, "y": 10.0, "z": 15.0}, - ] - } + ], + }, ] result_list = get_center_points_and_radius_for_each_ptm(ptm_list) ptm_1 = result_list[0] - + result_list = get_center_points_and_radius_for_each_ptm(ptm_list) ptm_1 = result_list[0] @@ -91,8 +85,8 @@ def test_get_center_points_and_radius_for_each_ptm(): assert pytest.approx(ptm_2["radius"], 1e-6) == 0.0 ptm_3 = result_list[2] - + expected_center = [10.0, 10.0, 10.0] np.testing.assert_allclose(ptm_3["center_point"], expected_center, atol=1e-6) - - assert pytest.approx(ptm_3["radius"], 1e-6) == 5.0 \ No newline at end of file + + assert pytest.approx(ptm_3["radius"], 1e-6) == 5.0 From d3e6085efc0991e11bb67a46d7d8c47bf00372fa Mon Sep 17 00:00:00 2001 From: 3dot141592 Date: Sun, 17 May 2026 08:52:54 +0900 Subject: [PATCH 12/22] finish vdw distance --- .../data_analysis/amino_acid_spheres.py | 8 ++-- .../data_analysis/geometry_operations.py | 37 +++++++++++++++---- 2 files changed, 33 insertions(+), 12 deletions(-) diff --git a/backend/protzilla/data_analysis/amino_acid_spheres.py b/backend/protzilla/data_analysis/amino_acid_spheres.py index ad7d00527..1cc4a48a3 100644 --- a/backend/protzilla/data_analysis/amino_acid_spheres.py +++ b/backend/protzilla/data_analysis/amino_acid_spheres.py @@ -1,13 +1,12 @@ from __future__ import annotations -import numpy as np import pandas as pd import trimesh from backend.protzilla.data_analysis.geometry_operations import ( calculate_center_point, extract_points_from_cif, - find_farthest_point, + find_farthest_point_vdw, mesh_to_polyhedron, ) @@ -29,14 +28,13 @@ def calculate_amino_acid_spheres( spheres = [] for residue_position in residue_positions: # TODO: This can certainly be made more efficient - residue_points = extract_points_from_cif( + residue_points, residue_elements = extract_points_from_cif( cif_df, residue_range=(residue_position, residue_position), chain_id=chain_id, ) center = calculate_center_point(residue_points) - furthest_point = find_farthest_point(residue_points, center) - radius = float(np.linalg.norm(furthest_point - center)) + _, radius = find_farthest_point_vdw(residue_points, residue_elements, center) sphere = trimesh.creation.icosphere(subdivisions=subdivisions, radius=radius) sphere.apply_translation(center) diff --git a/backend/protzilla/data_analysis/geometry_operations.py b/backend/protzilla/data_analysis/geometry_operations.py index 9605985e2..e29ed5d97 100644 --- a/backend/protzilla/data_analysis/geometry_operations.py +++ b/backend/protzilla/data_analysis/geometry_operations.py @@ -4,6 +4,8 @@ import pandas as pd import trimesh +from backend.protzilla.constants.van_der_waals import vdw_radii + COORDINATE_COLUMNS = [ "_atom_site.Cartn_x", "_atom_site.Cartn_y", @@ -27,7 +29,7 @@ def extract_points_from_cif( cif_df: pd.DataFrame, residue_range: tuple[int, int], chain_id: str | None = None, -) -> np.ndarray: +) -> tuple[np.ndarray, np.ndarray]: """ Extract Cartesian atom coordinates in a residue range from a CIF DataFrame. @@ -41,6 +43,7 @@ def extract_points_from_cif( required_columns = { "_atom_site.label_seq_id", + "_atom_site.type_symbol", "_atom_site.Cartn_x", "_atom_site.Cartn_y", "_atom_site.Cartn_z", @@ -70,8 +73,17 @@ def extract_points_from_cif( residue_ids = filtered_df["_atom_site.label_seq_id"].astype(int) filtered_df = filtered_df[(residue_ids >= start) & (residue_ids <= end)] + atom_data = filtered_df[ + [ + "_atom_site.Cartn_x", + "_atom_site.Cartn_y", + "_atom_site.Cartn_z", + "_atom_site.type_symbol", + ] + ].drop_duplicates() + points = ( - filtered_df[ + atom_data[ [ "_atom_site.Cartn_x", "_atom_site.Cartn_y", @@ -79,14 +91,20 @@ def extract_points_from_cif( ] ] .astype(float) - .drop_duplicates() + .to_numpy() + ) + elements = ( + atom_data["_atom_site.type_symbol"] + .astype(str) + .str.strip() + .str.capitalize() .to_numpy() ) if len(points) == 0: raise ValueError("No atom coordinates found.") - return points + return points, elements def build_convex_hull(points: np.ndarray) -> trimesh.Trimesh: @@ -146,7 +164,9 @@ def find_farthest_point(points: np.ndarray, reference_point: np.ndarray) -> np.n distances = np.linalg.norm(points - reference_point, axis=1) return points[np.argmax(distances)] -def find_farthest_point_fdw(points: np.ndarray, reference_point: np.ndarray) -> np.ndarray: +def find_farthest_point_vdw( + points: np.ndarray, elements: np.ndarray, reference_point: np.ndarray +) -> tuple[np.ndarray, float]: """ Calculate the point, whose fdw-"Bubble" is most distant to a reference point. """ @@ -159,8 +179,11 @@ def find_farthest_point_fdw(points: np.ndarray, reference_point: np.ndarray) -> f"Expected reference_point with shape (3,), got {reference_point.shape}." ) - distances = np.linalg.norm(points - reference_point, axis=1) - return points[np.argmax(distances)] + radii = np.array([vdw_radii[element] for element in elements], dtype=float) + + distances = np.linalg.norm(points - reference_point, axis=1) + radii + max_index = np.argmax(distances) + return points[max_index], float(distances[max_index]) def meshes_intersect( From a46eecec21a1db5a065cc29f7896ff9bd8260f70 Mon Sep 17 00:00:00 2001 From: 3dot141592 Date: Thu, 21 May 2026 12:52:14 +0900 Subject: [PATCH 13/22] ptm bubble implementation --- backend/main/views_helper.py | 6 ++- .../data_analysis/amino_acid_spheres.py | 37 +++++++++++++++++++ backend/protzilla/utilities/ptm_helper.py | 15 ++++++-- 3 files changed, 54 insertions(+), 4 deletions(-) diff --git a/backend/main/views_helper.py b/backend/main/views_helper.py index 934bcc3a4..2e13d4d1a 100644 --- a/backend/main/views_helper.py +++ b/backend/main/views_helper.py @@ -6,7 +6,10 @@ from typing import Optional, List, Dict from backend.protzilla.constants.paths import SETTINGS_PATH -from backend.protzilla.data_analysis.amino_acid_spheres import calculate_amino_acid_spheres +from backend.protzilla.data_analysis.amino_acid_spheres import ( + calculate_amino_acid_spheres, + calculate_ptm_spheres, +) from backend.protzilla.disk_operator import YamlOperator from backend.protzilla.steps import Step from backend.protzilla.step_manager import StepManager @@ -211,6 +214,7 @@ def create_visualization( result = {"structureEntryId": structure_entry_id, "cifString": cif_string} result["trimeshMeshes"] = calculate_amino_acid_spheres(cif_df) + result["trimeshMeshes"].extend(calculate_ptm_spheres(cif_df)) if crosslinking_df is not None: result["crosslinks"] = extract_relevant_crosslink_information(crosslinking_df) diff --git a/backend/protzilla/data_analysis/amino_acid_spheres.py b/backend/protzilla/data_analysis/amino_acid_spheres.py index 1cc4a48a3..f442a3045 100644 --- a/backend/protzilla/data_analysis/amino_acid_spheres.py +++ b/backend/protzilla/data_analysis/amino_acid_spheres.py @@ -1,5 +1,6 @@ from __future__ import annotations +import numpy as np import pandas as pd import trimesh @@ -9,6 +10,10 @@ find_farthest_point_vdw, mesh_to_polyhedron, ) +from backend.protzilla.utilities.ptm_helper import ( + get_all_ptm_atoms_with_coordinates, + get_center_points_and_radius_for_each_ptm, +) def calculate_amino_acid_spheres( @@ -49,3 +54,35 @@ def calculate_amino_acid_spheres( ) return spheres + + +def calculate_ptm_spheres( + cif_df: pd.DataFrame, + color: int = 0x00A6A6, + alpha: float = 0.35, + subdivisions: int = 1, +) -> list[dict]: + if "_chem_comp.mon_nstd_flag" not in cif_df.columns: + return [] + + ptms = get_all_ptm_atoms_with_coordinates(cif_df) + ptms = get_center_points_and_radius_for_each_ptm(ptms) + + spheres = [] + for ptm in ptms: + center = np.array(ptm["center_point"], dtype=float) + radius = float(ptm["radius"]) + + sphere = trimesh.creation.icosphere(subdivisions=subdivisions, radius=radius) + sphere.apply_translation(center) + + spheres.append( + { + "label": f"PTM {ptm['ptm_name']} {ptm['chain']}:{ptm['position']} sphere", + "mesh": mesh_to_polyhedron(sphere), + "color": color, + "alpha": alpha, + } + ) + + return spheres diff --git a/backend/protzilla/utilities/ptm_helper.py b/backend/protzilla/utilities/ptm_helper.py index 65d6a1229..b5b33b8c6 100644 --- a/backend/protzilla/utilities/ptm_helper.py +++ b/backend/protzilla/utilities/ptm_helper.py @@ -1,6 +1,8 @@ import pandas as pd import numpy as np +from backend.protzilla.data_analysis.geometry_operations import find_farthest_point_vdw + def get_all_ptm_atoms_with_coordinates(cif_df: pd.DataFrame) -> list: """ @@ -21,6 +23,7 @@ def get_all_ptm_atoms_with_coordinates(cif_df: pd.DataFrame) -> list: atoms_subset = atom_group[ [ "_atom_site.label_atom_id", + "_atom_site.type_symbol", "_atom_site.Cartn_x", "_atom_site.Cartn_y", "_atom_site.Cartn_z", @@ -28,11 +31,15 @@ def get_all_ptm_atoms_with_coordinates(cif_df: pd.DataFrame) -> list: ].rename( columns={ "_atom_site.label_atom_id": "atom_name", + "_atom_site.type_symbol": "element", "_atom_site.Cartn_x": "x", "_atom_site.Cartn_y": "y", "_atom_site.Cartn_z": "z", } ) + atoms_subset["element"] = ( + atoms_subset["element"].astype(str).str.strip().str.capitalize() + ) atoms_data = atoms_subset.to_dict(orient="records") @@ -82,13 +89,15 @@ def find_farthest_point(points: np.ndarray, reference_point: np.ndarray) -> np.n def get_center_points_and_radius_for_each_ptm(ptm_list: list) -> list: for ptm in ptm_list: coords = [[atom["x"], atom["y"], atom["z"]] for atom in ptm["atoms"]] + elements = [atom["element"] for atom in ptm["atoms"]] coords_array = np.array(coords, dtype=np.float32) + elements_array = np.array(elements) center_point = calculate_center_point(coords_array) - farthest_point = find_farthest_point(coords_array, center_point) - - radius = np.linalg.norm(farthest_point - center_point) + _, radius = find_farthest_point_vdw( + coords_array, elements_array, center_point + ) ptm["center_point"] = center_point.tolist() ptm["radius"] = float(radius) From bcf9a039a19a0c18ed167198c871462c83cf6028 Mon Sep 17 00:00:00 2001 From: 3dot141592 Date: Tue, 26 May 2026 13:24:57 +0900 Subject: [PATCH 14/22] PTM amino acid collision --- backend/main/views_helper.py | 4 +- .../data_analysis/amino_acid_spheres.py | 125 ++++++++++++++++-- .../data_analysis/geometry_operations.py | 23 +++- 3 files changed, 138 insertions(+), 14 deletions(-) diff --git a/backend/main/views_helper.py b/backend/main/views_helper.py index 2e13d4d1a..f9d648d37 100644 --- a/backend/main/views_helper.py +++ b/backend/main/views_helper.py @@ -213,7 +213,9 @@ def create_visualization( result = {"structureEntryId": structure_entry_id, "cifString": cif_string} - result["trimeshMeshes"] = calculate_amino_acid_spheres(cif_df) + result["trimeshMeshes"] = calculate_amino_acid_spheres( + cif_df, only_intersecting_ptms=True, ignored_neighbors=0 + ) result["trimeshMeshes"].extend(calculate_ptm_spheres(cif_df)) if crosslinking_df is not None: diff --git a/backend/protzilla/data_analysis/amino_acid_spheres.py b/backend/protzilla/data_analysis/amino_acid_spheres.py index f442a3045..41fd8add9 100644 --- a/backend/protzilla/data_analysis/amino_acid_spheres.py +++ b/backend/protzilla/data_analysis/amino_acid_spheres.py @@ -8,7 +8,9 @@ calculate_center_point, extract_points_from_cif, find_farthest_point_vdw, + find_intersecting_spheres, mesh_to_polyhedron, + resolve_chain_column, ) from backend.protzilla.utilities.ptm_helper import ( get_all_ptm_atoms_with_coordinates, @@ -16,15 +18,16 @@ ) -def calculate_amino_acid_spheres( - cif_df: pd.DataFrame, - chain_id: str | None = None, - color: int = 0xFF8C00, - alpha: float = 0.25, - subdivisions: int = 1, +def _calculate_residue_spheres( + cif_df: pd.DataFrame, chain_id: str | None = None ) -> list[dict]: + residue_df = cif_df + chain_column = resolve_chain_column(cif_df) + if chain_id is not None and chain_column is not None: + residue_df = cif_df[cif_df[chain_column] == chain_id] + residue_positions = ( - pd.to_numeric(cif_df["_atom_site.label_seq_id"], errors="coerce") + pd.to_numeric(residue_df["_atom_site.label_seq_id"], errors="coerce") .dropna() .astype(int) .drop_duplicates() @@ -32,7 +35,7 @@ def calculate_amino_acid_spheres( ) spheres = [] - for residue_position in residue_positions: # TODO: This can certainly be made more efficient + for residue_position in residue_positions: residue_points, residue_elements = extract_points_from_cif( cif_df, residue_range=(residue_position, residue_position), @@ -41,12 +44,51 @@ def calculate_amino_acid_spheres( center = calculate_center_point(residue_points) _, radius = find_farthest_point_vdw(residue_points, residue_elements, center) - sphere = trimesh.creation.icosphere(subdivisions=subdivisions, radius=radius) - sphere.apply_translation(center) + spheres.append( + { + "chain": chain_id, + "position": residue_position, + "center": center.tolist(), + "radius": float(radius), + } + ) + + return spheres + + +def calculate_amino_acid_spheres( + cif_df: pd.DataFrame, + chain_id: str | None = None, + color: int = 0xFF8C00, + alpha: float = 0.25, + subdivisions: int = 1, + only_intersecting_ptms: bool = False, + ignored_neighbors: int = 0, +) -> list[dict]: + if only_intersecting_ptms: + residue_spheres = [] + seen_residues = set() + for collision in find_ptm_amino_acid_sphere_collisions( + cif_df, ignored_neighbors=ignored_neighbors, chain_id=chain_id + ): + for residue_sphere in collision["collisions"]: + residue_key = (residue_sphere["chain"], residue_sphere["position"]) + if residue_key not in seen_residues: + seen_residues.add(residue_key) + residue_spheres.append(residue_sphere) + else: + residue_spheres = _calculate_residue_spheres(cif_df, chain_id) + + spheres = [] + for residue_sphere in residue_spheres: + sphere = trimesh.creation.icosphere( + subdivisions=subdivisions, radius=residue_sphere["radius"] + ) + sphere.apply_translation(residue_sphere["center"]) spheres.append( { - "label": f"Residue {residue_position} sphere", + "label": f"Residue {residue_sphere['position']} sphere", "mesh": mesh_to_polyhedron(sphere), "color": color, "alpha": alpha, @@ -86,3 +128,64 @@ def calculate_ptm_spheres( ) return spheres + + +def find_ptm_amino_acid_sphere_collisions( + cif_df: pd.DataFrame, + ignored_neighbors: int = 0, + chain_id: str | None = None, +) -> list[dict]: + if "_chem_comp.mon_nstd_flag" not in cif_df.columns: + return [] + + amino_acid_df = cif_df[cif_df["_chem_comp.mon_nstd_flag"] == True] + chain_column = resolve_chain_column(amino_acid_df) + if chain_id is not None: + chains = [chain_id] + elif chain_column is not None: + chains = amino_acid_df[chain_column].dropna().drop_duplicates() + else: + chains = [None] + amino_acid_spheres = [ + sphere + for chain in chains + for sphere in _calculate_residue_spheres(amino_acid_df, chain) + ] + + ptms = get_center_points_and_radius_for_each_ptm( + get_all_ptm_atoms_with_coordinates(cif_df) + ) + if chain_id is not None: + ptms = [ptm for ptm in ptms if ptm["chain"] == chain_id] + + collisions = [] + for ptm in ptms: + ptm_sphere = { + "ptm_name": ptm["ptm_name"], + "chain": ptm["chain"], + "position": int(ptm["position"]), + "center": ptm["center_point"], + "radius": float(ptm["radius"]), + } + checked_spheres = [ + sphere + for sphere in amino_acid_spheres + if not ( + sphere["chain"] == ptm_sphere["chain"] + and 0 + < abs(sphere["position"] - ptm_sphere["position"]) + <= ignored_neighbors + ) + ] + intersecting_spheres = find_intersecting_spheres( + checked_spheres, ptm_sphere + ) + + collisions.append( + { + **ptm_sphere, + "collisions": intersecting_spheres, + } + ) + + return collisions diff --git a/backend/protzilla/data_analysis/geometry_operations.py b/backend/protzilla/data_analysis/geometry_operations.py index e29ed5d97..58e028d0e 100644 --- a/backend/protzilla/data_analysis/geometry_operations.py +++ b/backend/protzilla/data_analysis/geometry_operations.py @@ -13,7 +13,7 @@ ] -def _resolve_chain_column(cif_df: pd.DataFrame) -> str | None: +def resolve_chain_column(cif_df: pd.DataFrame) -> str | None: """ Return the preferred chain identifier column if present in the CIF DataFrame. """ @@ -62,7 +62,7 @@ def extract_points_from_cif( filtered_df = cif_df.copy() - chain_column = _resolve_chain_column(filtered_df) + chain_column = resolve_chain_column(filtered_df) if chain_id is not None: if chain_column is None: raise ValueError( @@ -186,6 +186,25 @@ def find_farthest_point_vdw( return points[max_index], float(distances[max_index]) +def find_intersecting_spheres( + spheres: list[dict], + reference_sphere: dict, +) -> list[dict]: + reference_center = np.array(reference_sphere["center"], dtype=float) + reference_radius = float(reference_sphere["radius"]) + + intersecting_spheres = [] + for sphere in spheres: + center = np.array(sphere["center"], dtype=float) + radius = float(sphere["radius"]) + distance = np.linalg.norm(center - reference_center) + + if distance <= reference_radius + radius: + intersecting_spheres.append(sphere) + + return intersecting_spheres + + def meshes_intersect( mesh_a: trimesh.Trimesh, mesh_b: trimesh.Trimesh, From fdb5a709410413e882a3060c608d0c365b8e71f3 Mon Sep 17 00:00:00 2001 From: 3dot141592 Date: Tue, 26 May 2026 16:04:08 +0900 Subject: [PATCH 15/22] black --- backend/main/views_helper.py | 4 ++-- backend/protzilla/constants/van_der_waals.py | 3 +-- .../protzilla/data_analysis/amino_acid_spheres.py | 4 +--- .../protzilla/data_analysis/geometry_operations.py | 13 +++++++++---- backend/protzilla/utilities/ptm_helper.py | 4 +--- 5 files changed, 14 insertions(+), 14 deletions(-) diff --git a/backend/main/views_helper.py b/backend/main/views_helper.py index f9d648d37..d686557b9 100644 --- a/backend/main/views_helper.py +++ b/backend/main/views_helper.py @@ -212,12 +212,12 @@ def create_visualization( cif_string = "" result = {"structureEntryId": structure_entry_id, "cifString": cif_string} - + result["trimeshMeshes"] = calculate_amino_acid_spheres( cif_df, only_intersecting_ptms=True, ignored_neighbors=0 ) result["trimeshMeshes"].extend(calculate_ptm_spheres(cif_df)) - + if crosslinking_df is not None: result["crosslinks"] = extract_relevant_crosslink_information(crosslinking_df) diff --git a/backend/protzilla/constants/van_der_waals.py b/backend/protzilla/constants/van_der_waals.py index a7e25ad8c..eb3519e6b 100644 --- a/backend/protzilla/constants/van_der_waals.py +++ b/backend/protzilla/constants/van_der_waals.py @@ -2,7 +2,6 @@ """ Van der Waals radius in Ångström """ - "H": 1.20, "He": 1.43, "Li": 2.12, @@ -121,4 +120,4 @@ "Lv": None, "Ts": None, "Og": None, -} \ No newline at end of file +} diff --git a/backend/protzilla/data_analysis/amino_acid_spheres.py b/backend/protzilla/data_analysis/amino_acid_spheres.py index 41fd8add9..5f5731313 100644 --- a/backend/protzilla/data_analysis/amino_acid_spheres.py +++ b/backend/protzilla/data_analysis/amino_acid_spheres.py @@ -177,9 +177,7 @@ def find_ptm_amino_acid_sphere_collisions( <= ignored_neighbors ) ] - intersecting_spheres = find_intersecting_spheres( - checked_spheres, ptm_sphere - ) + intersecting_spheres = find_intersecting_spheres(checked_spheres, ptm_sphere) collisions.append( { diff --git a/backend/protzilla/data_analysis/geometry_operations.py b/backend/protzilla/data_analysis/geometry_operations.py index 58e028d0e..b3c677edd 100644 --- a/backend/protzilla/data_analysis/geometry_operations.py +++ b/backend/protzilla/data_analysis/geometry_operations.py @@ -138,7 +138,7 @@ def point_cloud_to_polyhedron(points: np.ndarray) -> dict: def calculate_center_point(points: np.ndarray) -> np.ndarray: """ - Calculate the center point of a point cloud. + Calculate the center point of a point cloud. """ if points.ndim != 2 or points.shape[1] != 3: raise ValueError(f"Expected points with shape (n, 3), got {points.shape}.") @@ -150,12 +150,14 @@ def calculate_center_point(points: np.ndarray) -> np.ndarray: def find_farthest_point(points: np.ndarray, reference_point: np.ndarray) -> np.ndarray: """ - Calculate the farthest point inside a point cloud from a reference point. + Calculate the farthest point inside a point cloud from a reference point. """ if points.ndim != 2 or points.shape[1] != 3: raise ValueError(f"Expected points with shape (n, 3), got {points.shape}.") if len(points) == 0: - raise ValueError("At least one point is required to calculate a maximum distance.") + raise ValueError( + "At least one point is required to calculate a maximum distance." + ) if reference_point.shape != (3,): raise ValueError( f"Expected reference_point with shape (3,), got {reference_point.shape}." @@ -164,6 +166,7 @@ def find_farthest_point(points: np.ndarray, reference_point: np.ndarray) -> np.n distances = np.linalg.norm(points - reference_point, axis=1) return points[np.argmax(distances)] + def find_farthest_point_vdw( points: np.ndarray, elements: np.ndarray, reference_point: np.ndarray ) -> tuple[np.ndarray, float]: @@ -173,7 +176,9 @@ def find_farthest_point_vdw( if points.ndim != 2 or points.shape[1] != 3: raise ValueError(f"Expected points with shape (n, 3), got {points.shape}.") if len(points) == 0: - raise ValueError("At least one point is required to calculate a maximum distance.") + raise ValueError( + "At least one point is required to calculate a maximum distance." + ) if reference_point.shape != (3,): raise ValueError( f"Expected reference_point with shape (3,), got {reference_point.shape}." diff --git a/backend/protzilla/utilities/ptm_helper.py b/backend/protzilla/utilities/ptm_helper.py index b5b33b8c6..e42c5f103 100644 --- a/backend/protzilla/utilities/ptm_helper.py +++ b/backend/protzilla/utilities/ptm_helper.py @@ -95,9 +95,7 @@ def get_center_points_and_radius_for_each_ptm(ptm_list: list) -> list: elements_array = np.array(elements) center_point = calculate_center_point(coords_array) - _, radius = find_farthest_point_vdw( - coords_array, elements_array, center_point - ) + _, radius = find_farthest_point_vdw(coords_array, elements_array, center_point) ptm["center_point"] = center_point.tolist() ptm["radius"] = float(radius) From ba08b5243e0ab208200c539d78a580c49925350d Mon Sep 17 00:00:00 2001 From: 3dot141592 Date: Thu, 28 May 2026 10:52:12 +0900 Subject: [PATCH 16/22] Add PTM collision output to validation step --- backend/main/views_helper.py | 13 +++-- .../protzilla/data_analysis/ptm_analysis.py | 49 ++++++++++++++++--- backend/protzilla/methods/data_analysis.py | 11 ++++- 3 files changed, 61 insertions(+), 12 deletions(-) diff --git a/backend/main/views_helper.py b/backend/main/views_helper.py index d686557b9..a45959077 100644 --- a/backend/main/views_helper.py +++ b/backend/main/views_helper.py @@ -191,6 +191,8 @@ def create_visualization( cif_df: pd.DataFrame, structure_entry_id: str, crosslinking_df: Optional[pd.DataFrame] = None, + include_ptm_spheres: bool = False, + ignored_neighbors: int = 0, ) -> dict: """ Create visualization data, by packaging a mmCIF string (converted from a CIF DataFrame) with its structure entry ID. @@ -213,10 +215,13 @@ def create_visualization( result = {"structureEntryId": structure_entry_id, "cifString": cif_string} - result["trimeshMeshes"] = calculate_amino_acid_spheres( - cif_df, only_intersecting_ptms=True, ignored_neighbors=0 - ) - result["trimeshMeshes"].extend(calculate_ptm_spheres(cif_df)) + if include_ptm_spheres: + result["trimeshMeshes"] = calculate_amino_acid_spheres( + cif_df, + only_intersecting_ptms=True, + ignored_neighbors=int(ignored_neighbors), + ) + result["trimeshMeshes"].extend(calculate_ptm_spheres(cif_df)) if crosslinking_df is not None: result["crosslinks"] = extract_relevant_crosslink_information(crosslinking_df) diff --git a/backend/protzilla/data_analysis/ptm_analysis.py b/backend/protzilla/data_analysis/ptm_analysis.py index 8665eaaea..aa07085ae 100644 --- a/backend/protzilla/data_analysis/ptm_analysis.py +++ b/backend/protzilla/data_analysis/ptm_analysis.py @@ -3,10 +3,10 @@ import re from backend.protzilla.utilities.transform_dfs import long_to_wide -from backend.protzilla.utilities.ptm_helper import ( - get_all_ptm_atoms_with_coordinates, - get_center_points_and_radius_for_each_ptm, +from backend.protzilla.data_analysis.amino_acid_spheres import ( + find_ptm_amino_acid_sphere_collisions, ) +from backend.protzilla.steps import OutputItem, OutputType def ptms_per_sample(psm_df: pd.DataFrame) -> dict: @@ -97,7 +97,42 @@ def from_string(mod_string: str) -> tuple[int, str]: return amount, name -def ptm_validation(cif_df: pd.DataFrame): - ptm_list = get_all_ptm_atoms_with_coordinates(cif_df) - ptm_list = get_center_points_and_radius_for_each_ptm(ptm_list) - return {"ptm_validation_df": pd.DataFrame(ptm_list)} +def ptm_validation( + cif_df: pd.DataFrame, + structure_metadata_df: pd.DataFrame, + ignored_neighbors: int = 0, +): + ignored_neighbors = int(ignored_neighbors) + collisions = find_ptm_amino_acid_sphere_collisions( + cif_df, ignored_neighbors=ignored_neighbors + ) + collision_columns = [ + "ptm_name", + "ptm_chain", + "ptm_position", + "amino_acid_chain", + "amino_acid_position", + ] + collision_rows = [ + { + "ptm_name": collision["ptm_name"], + "ptm_chain": collision["chain"], + "ptm_position": collision["position"], + "amino_acid_chain": amino_acid["chain"], + "amino_acid_position": amino_acid["position"], + } + for collision in collisions + for amino_acid in collision["collisions"] + ] + data_for_visualization = { + "structure_entry_id": structure_metadata_df["entry_id"].iloc[0], + "cif_df": cif_df, + "include_ptm_spheres": True, + "ignored_neighbors": ignored_neighbors, + } + return { + "ptm_validation_df": pd.DataFrame(collision_rows, columns=collision_columns), + "visualization": OutputItem( + output_type=OutputType.VISUALIZATION, value=data_for_visualization + ), + } diff --git a/backend/protzilla/methods/data_analysis.py b/backend/protzilla/methods/data_analysis.py index a494d159b..19e06b2bb 100644 --- a/backend/protzilla/methods/data_analysis.py +++ b/backend/protzilla/methods/data_analysis.py @@ -2488,5 +2488,14 @@ class PtmValidation(DataAnalysisStep): def create_form(self): return Form( label="PTM Validation", - input_fields=[], + input_fields=[ + NumberField( + name="ignored_neighbors", + label="Number of PTM residue neighbors to ignore", + value=0, + min=0, + step=1, + hasStepButtons=True, + ) + ], ) From 3df865f2d47326c100cd1d7b62a044461aa964c5 Mon Sep 17 00:00:00 2001 From: 3dot141592 Date: Mon, 15 Jun 2026 10:16:45 +0200 Subject: [PATCH 17/22] PR: Use CIF column constants --- backend/protzilla/constants/cif_columns.py | 1 + .../data_analysis/amino_acid_spheres.py | 12 +++-- .../data_analysis/crosslinking_validation.py | 19 ++++---- .../data_analysis/geometry_operations.py | 46 ++++++------------- backend/protzilla/utilities/ptm_helper.py | 30 ++++++------ 5 files changed, 51 insertions(+), 57 deletions(-) diff --git a/backend/protzilla/constants/cif_columns.py b/backend/protzilla/constants/cif_columns.py index 5715f2f32..2d0f65afd 100644 --- a/backend/protzilla/constants/cif_columns.py +++ b/backend/protzilla/constants/cif_columns.py @@ -30,6 +30,7 @@ class ATOM_SITE_COLUMNS(StrEnum): AUTH_ASYM_ID = f"{ATOM_SITE_PREFIX}auth_asym_id" AUTH_ATOM_ID = f"{ATOM_SITE_PREFIX}auth_atom_id" PDBX_PDB_MODEL_NUM = f"{ATOM_SITE_PREFIX}pdbx_PDB_model_num" + PDBX_SIFTS_XREF_DB_ACC = f"{ATOM_SITE_PREFIX}pdbx_sifts_xref_db_acc" ATOM_SITE_LABEL_COMP_ID = ATOM_SITE_COLUMNS.LABEL_COMP_ID diff --git a/backend/protzilla/data_analysis/amino_acid_spheres.py b/backend/protzilla/data_analysis/amino_acid_spheres.py index 5f5731313..ff3150868 100644 --- a/backend/protzilla/data_analysis/amino_acid_spheres.py +++ b/backend/protzilla/data_analysis/amino_acid_spheres.py @@ -4,6 +4,10 @@ import pandas as pd import trimesh +from backend.protzilla.constants.cif_columns import ( + ATOM_SITE_COLUMNS, + CHEM_COMP_COLUMNS, +) from backend.protzilla.data_analysis.geometry_operations import ( calculate_center_point, extract_points_from_cif, @@ -27,7 +31,7 @@ def _calculate_residue_spheres( residue_df = cif_df[cif_df[chain_column] == chain_id] residue_positions = ( - pd.to_numeric(residue_df["_atom_site.label_seq_id"], errors="coerce") + pd.to_numeric(residue_df[ATOM_SITE_COLUMNS.LABEL_SEQ_ID], errors="coerce") .dropna() .astype(int) .drop_duplicates() @@ -104,7 +108,7 @@ def calculate_ptm_spheres( alpha: float = 0.35, subdivisions: int = 1, ) -> list[dict]: - if "_chem_comp.mon_nstd_flag" not in cif_df.columns: + if CHEM_COMP_COLUMNS.MON_NSTD_FLAG not in cif_df.columns: return [] ptms = get_all_ptm_atoms_with_coordinates(cif_df) @@ -135,10 +139,10 @@ def find_ptm_amino_acid_sphere_collisions( ignored_neighbors: int = 0, chain_id: str | None = None, ) -> list[dict]: - if "_chem_comp.mon_nstd_flag" not in cif_df.columns: + if CHEM_COMP_COLUMNS.MON_NSTD_FLAG not in cif_df.columns: return [] - amino_acid_df = cif_df[cif_df["_chem_comp.mon_nstd_flag"] == True] + amino_acid_df = cif_df[cif_df[CHEM_COMP_COLUMNS.MON_NSTD_FLAG] == True] chain_column = resolve_chain_column(amino_acid_df) if chain_id is not None: chains = [chain_id] diff --git a/backend/protzilla/data_analysis/crosslinking_validation.py b/backend/protzilla/data_analysis/crosslinking_validation.py index 4f02ee3bb..e9deb5f57 100644 --- a/backend/protzilla/data_analysis/crosslinking_validation.py +++ b/backend/protzilla/data_analysis/crosslinking_validation.py @@ -34,6 +34,7 @@ PLOT_PRIMARY_COLOR, PLOT_SECONDARY_COLOR, ) +from backend.protzilla.constants.cif_columns import ATOM_SITE_COLUMNS def get_reactive_atom_of_amino_acid_residue(amino_acid_type: str) -> str: @@ -70,14 +71,14 @@ def get_coordinates_of_atom_crosslinker_bound_to( """ relevant_atom = get_reactive_atom_of_amino_acid_residue(amino_acid_type) - seq_ids = pd.to_numeric(cif_df["_atom_site.label_seq_id"], errors="coerce") + seq_ids = pd.to_numeric(cif_df[ATOM_SITE_COLUMNS.LABEL_SEQ_ID], errors="coerce") # Filter to the exact reactive atom of the amino acid residue # where the crosslinker is bound (e.g. CA at position 45) cif_df = cif_df[ - (cif_df["_atom_site.label_atom_id"] == relevant_atom) + (cif_df[ATOM_SITE_COLUMNS.LABEL_ATOM_ID] == relevant_atom) & (seq_ids == amino_acid_position_where_crosslinker_bound) - & (cif_df["_atom_site.auth_asym_id"] == chain_id) + & (cif_df[ATOM_SITE_COLUMNS.AUTH_ASYM_ID] == chain_id) ] if cif_df.empty: @@ -87,9 +88,9 @@ def get_coordinates_of_atom_crosslinker_bound_to( row = cif_df.iloc[0] - x = float(row["_atom_site.Cartn_x"]) - y = float(row["_atom_site.Cartn_y"]) - z = float(row["_atom_site.Cartn_z"]) + x = float(row[ATOM_SITE_COLUMNS.CARTN_X]) + y = float(row[ATOM_SITE_COLUMNS.CARTN_Y]) + z = float(row[ATOM_SITE_COLUMNS.CARTN_Z]) return x, y, z @@ -308,7 +309,7 @@ def get_chains( return [] target_ids_as_strings = [str(i) for i in target_ids] relevant_df = cif_df[cif_df[id_column_name].astype(str).isin(target_ids_as_strings)] - chain_ids = relevant_df["_atom_site.auth_asym_id"].dropna().unique().tolist() + chain_ids = relevant_df[ATOM_SITE_COLUMNS.AUTH_ASYM_ID].dropna().unique().tolist() return chain_ids @@ -393,7 +394,7 @@ def monomer_validation( cif_df=cif_df, amino_acid_sequences_df=amino_acid_sequences_df, valid_ids=valid_ids, - id_column_name="_atom_site.pdbx_sifts_xref_db_acc", + id_column_name=ATOM_SITE_COLUMNS.PDBX_SIFTS_XREF_DB_ACC, structures_to_validate=[protein_id], ) @@ -491,7 +492,7 @@ def multimer_validation( cif_df=cif_df, amino_acid_sequences_df=amino_acid_sequences_df, valid_ids=valid_ids, - id_column_name="_atom_site.label_entity_id", + id_column_name=ATOM_SITE_COLUMNS.LABEL_ENTITY_ID, structures_to_validate=structures_to_validate, ) diff --git a/backend/protzilla/data_analysis/geometry_operations.py b/backend/protzilla/data_analysis/geometry_operations.py index b3c677edd..1d7428d1f 100644 --- a/backend/protzilla/data_analysis/geometry_operations.py +++ b/backend/protzilla/data_analysis/geometry_operations.py @@ -4,12 +4,13 @@ import pandas as pd import trimesh +from backend.protzilla.constants.cif_columns import ATOM_SITE_COLUMNS from backend.protzilla.constants.van_der_waals import vdw_radii COORDINATE_COLUMNS = [ - "_atom_site.Cartn_x", - "_atom_site.Cartn_y", - "_atom_site.Cartn_z", + ATOM_SITE_COLUMNS.CARTN_X, + ATOM_SITE_COLUMNS.CARTN_Y, + ATOM_SITE_COLUMNS.CARTN_Z, ] @@ -18,10 +19,10 @@ def resolve_chain_column(cif_df: pd.DataFrame) -> str | None: Return the preferred chain identifier column if present in the CIF DataFrame. """ - if "_atom_site.label_asym_id" in cif_df.columns: - return "_atom_site.label_asym_id" - if "_atom_site.auth_asym_id" in cif_df.columns: - return "_atom_site.auth_asym_id" + if ATOM_SITE_COLUMNS.LABEL_ASYM_ID in cif_df.columns: + return ATOM_SITE_COLUMNS.LABEL_ASYM_ID + if ATOM_SITE_COLUMNS.AUTH_ASYM_ID in cif_df.columns: + return ATOM_SITE_COLUMNS.AUTH_ASYM_ID return None @@ -42,11 +43,9 @@ def extract_points_from_cif( """ required_columns = { - "_atom_site.label_seq_id", - "_atom_site.type_symbol", - "_atom_site.Cartn_x", - "_atom_site.Cartn_y", - "_atom_site.Cartn_z", + ATOM_SITE_COLUMNS.LABEL_SEQ_ID, + ATOM_SITE_COLUMNS.TYPE_SYMBOL, + *COORDINATE_COLUMNS, } missing_columns = sorted(required_columns - set(cif_df.columns)) if missing_columns: @@ -70,31 +69,16 @@ def extract_points_from_cif( ) filtered_df = filtered_df[filtered_df[chain_column] == chain_id] - residue_ids = filtered_df["_atom_site.label_seq_id"].astype(int) + residue_ids = filtered_df[ATOM_SITE_COLUMNS.LABEL_SEQ_ID].astype(int) filtered_df = filtered_df[(residue_ids >= start) & (residue_ids <= end)] atom_data = filtered_df[ - [ - "_atom_site.Cartn_x", - "_atom_site.Cartn_y", - "_atom_site.Cartn_z", - "_atom_site.type_symbol", - ] + [*COORDINATE_COLUMNS, ATOM_SITE_COLUMNS.TYPE_SYMBOL] ].drop_duplicates() - points = ( - atom_data[ - [ - "_atom_site.Cartn_x", - "_atom_site.Cartn_y", - "_atom_site.Cartn_z", - ] - ] - .astype(float) - .to_numpy() - ) + points = atom_data[COORDINATE_COLUMNS].astype(float).to_numpy() elements = ( - atom_data["_atom_site.type_symbol"] + atom_data[ATOM_SITE_COLUMNS.TYPE_SYMBOL] .astype(str) .str.strip() .str.capitalize() diff --git a/backend/protzilla/utilities/ptm_helper.py b/backend/protzilla/utilities/ptm_helper.py index e42c5f103..be36a8738 100644 --- a/backend/protzilla/utilities/ptm_helper.py +++ b/backend/protzilla/utilities/ptm_helper.py @@ -1,6 +1,10 @@ import pandas as pd import numpy as np +from backend.protzilla.constants.cif_columns import ( + ATOM_SITE_COLUMNS, + CHEM_COMP_COLUMNS, +) from backend.protzilla.data_analysis.geometry_operations import find_farthest_point_vdw @@ -11,30 +15,30 @@ def get_all_ptm_atoms_with_coordinates(cif_df: pd.DataFrame) -> list: """ ptm_list = [] - is_non_standard = cif_df["_chem_comp.mon_nstd_flag"] == False + is_non_standard = cif_df[CHEM_COMP_COLUMNS.MON_NSTD_FLAG] == False ptm_atoms_df = cif_df[is_non_standard] grouped_ptms = ptm_atoms_df.groupby( - ["_atom_site.label_asym_id", "_atom_site.label_seq_id"] + [ATOM_SITE_COLUMNS.LABEL_ASYM_ID, ATOM_SITE_COLUMNS.LABEL_SEQ_ID] ) for (chain_id, seq_id), atom_group in grouped_ptms: - comp_id = atom_group["_atom_site.label_comp_id"].iloc[0] + comp_id = atom_group[ATOM_SITE_COLUMNS.LABEL_COMP_ID].iloc[0] atoms_subset = atom_group[ [ - "_atom_site.label_atom_id", - "_atom_site.type_symbol", - "_atom_site.Cartn_x", - "_atom_site.Cartn_y", - "_atom_site.Cartn_z", + ATOM_SITE_COLUMNS.LABEL_ATOM_ID, + ATOM_SITE_COLUMNS.TYPE_SYMBOL, + ATOM_SITE_COLUMNS.CARTN_X, + ATOM_SITE_COLUMNS.CARTN_Y, + ATOM_SITE_COLUMNS.CARTN_Z, ] ].rename( columns={ - "_atom_site.label_atom_id": "atom_name", - "_atom_site.type_symbol": "element", - "_atom_site.Cartn_x": "x", - "_atom_site.Cartn_y": "y", - "_atom_site.Cartn_z": "z", + ATOM_SITE_COLUMNS.LABEL_ATOM_ID: "atom_name", + ATOM_SITE_COLUMNS.TYPE_SYMBOL: "element", + ATOM_SITE_COLUMNS.CARTN_X: "x", + ATOM_SITE_COLUMNS.CARTN_Y: "y", + ATOM_SITE_COLUMNS.CARTN_Z: "z", } ) atoms_subset["element"] = ( From 703bc294958fbbdae8bacc2bafa9600e0323d294 Mon Sep 17 00:00:00 2001 From: 3dot141592 <128993914+3dot141592@users.noreply.github.com> Date: Mon, 15 Jun 2026 10:19:39 +0200 Subject: [PATCH 18/22] PR: Code style Co-authored-by: Tarek --- backend/protzilla/data_analysis/amino_acid_spheres.py | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/backend/protzilla/data_analysis/amino_acid_spheres.py b/backend/protzilla/data_analysis/amino_acid_spheres.py index ff3150868..213d2bb43 100644 --- a/backend/protzilla/data_analysis/amino_acid_spheres.py +++ b/backend/protzilla/data_analysis/amino_acid_spheres.py @@ -111,8 +111,9 @@ def calculate_ptm_spheres( if CHEM_COMP_COLUMNS.MON_NSTD_FLAG not in cif_df.columns: return [] - ptms = get_all_ptm_atoms_with_coordinates(cif_df) - ptms = get_center_points_and_radius_for_each_ptm(ptms) + ptms = get_center_points_and_radius_for_each_ptm( + get_all_ptm_atoms_with_coordinates(cif_df) + ) spheres = [] for ptm in ptms: From f174295fe47a84b5b4ef022f93dc3497fe2d974c Mon Sep 17 00:00:00 2001 From: 3dot141592 Date: Mon, 15 Jun 2026 10:34:29 +0200 Subject: [PATCH 19/22] PR: rename ptm_validation_df to ptm_collisions_df --- backend/protzilla/data_analysis/ptm_analysis.py | 2 +- backend/protzilla/methods/data_analysis.py | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/backend/protzilla/data_analysis/ptm_analysis.py b/backend/protzilla/data_analysis/ptm_analysis.py index aa07085ae..de2789c2d 100644 --- a/backend/protzilla/data_analysis/ptm_analysis.py +++ b/backend/protzilla/data_analysis/ptm_analysis.py @@ -131,7 +131,7 @@ def ptm_validation( "ignored_neighbors": ignored_neighbors, } return { - "ptm_validation_df": pd.DataFrame(collision_rows, columns=collision_columns), + "ptm_collisions_df": pd.DataFrame(collision_rows, columns=collision_columns), "visualization": OutputItem( output_type=OutputType.VISUALIZATION, value=data_for_visualization ), diff --git a/backend/protzilla/methods/data_analysis.py b/backend/protzilla/methods/data_analysis.py index 19e06b2bb..671a2b72b 100644 --- a/backend/protzilla/methods/data_analysis.py +++ b/backend/protzilla/methods/data_analysis.py @@ -2482,7 +2482,7 @@ class PtmValidation(DataAnalysisStep): display_name = "PTM Validation" operation = "Peptide analysis" method_description = "Validates PTMs in protein structure predictions." - output_keys = ["ptm_validation_df"] + output_keys = ["ptm_collisions_df"] calc_method = staticmethod(ptm_validation) def create_form(self): From ae2e5e0c38fa0a5295a639740c44c0a1044aa8a8 Mon Sep 17 00:00:00 2001 From: 3dot141592 Date: Mon, 15 Jun 2026 10:43:36 +0200 Subject: [PATCH 20/22] PR: warning when no PTMs found --- backend/protzilla/data_analysis/ptm_analysis.py | 10 ++++++++++ 1 file changed, 10 insertions(+) diff --git a/backend/protzilla/data_analysis/ptm_analysis.py b/backend/protzilla/data_analysis/ptm_analysis.py index de2789c2d..983332899 100644 --- a/backend/protzilla/data_analysis/ptm_analysis.py +++ b/backend/protzilla/data_analysis/ptm_analysis.py @@ -1,6 +1,7 @@ import numpy as np import pandas as pd import re +import logging from backend.protzilla.utilities.transform_dfs import long_to_wide from backend.protzilla.data_analysis.amino_acid_spheres import ( @@ -106,6 +107,14 @@ def ptm_validation( collisions = find_ptm_amino_acid_sphere_collisions( cif_df, ignored_neighbors=ignored_neighbors ) + messages = [] + if not collisions: + messages.append( + dict( + level=logging.WARNING, + msg="No PTMs were found in the provided structure or PTM annotations are missing.", + ) + ) collision_columns = [ "ptm_name", "ptm_chain", @@ -135,4 +144,5 @@ def ptm_validation( "visualization": OutputItem( output_type=OutputType.VISUALIZATION, value=data_for_visualization ), + "messages": messages, } From bf09aaa326e65eead11b25d51cf5a49a4b9fb0d6 Mon Sep 17 00:00:00 2001 From: 3dot141592 Date: Mon, 15 Jun 2026 12:02:58 +0200 Subject: [PATCH 21/22] Add docstrings --- .../data_analysis/geometry_operations.py | 40 +++++++++++++++++++ 1 file changed, 40 insertions(+) diff --git a/backend/protzilla/data_analysis/geometry_operations.py b/backend/protzilla/data_analysis/geometry_operations.py index 1d7428d1f..1ef9bb5d8 100644 --- a/backend/protzilla/data_analysis/geometry_operations.py +++ b/backend/protzilla/data_analysis/geometry_operations.py @@ -17,6 +17,9 @@ def resolve_chain_column(cif_df: pd.DataFrame) -> str | None: """ Return the preferred chain identifier column if present in the CIF DataFrame. + + :param cif_df: DataFrame containing mmCIF atom_site data + :return: label chain column, author chain column as fallback, or None """ if ATOM_SITE_COLUMNS.LABEL_ASYM_ID in cif_df.columns: @@ -94,6 +97,11 @@ def extract_points_from_cif( def build_convex_hull(points: np.ndarray) -> trimesh.Trimesh: """ Build a 3D convex hull mesh from a point cloud. + + :param points: point coordinates with shape (n, 3) + :return: convex hull + :raises ValueError: if the param points have an invalid shape or contain fewer than + four points """ if points.ndim != 2 or points.shape[1] != 3: @@ -123,6 +131,10 @@ def point_cloud_to_polyhedron(points: np.ndarray) -> dict: def calculate_center_point(points: np.ndarray) -> np.ndarray: """ Calculate the center point of a point cloud. + + :param points: point coordinates with shape (n, 3) + :return: arithmetic mean of all points + :raises ValueError: if the param points have an invalid shape or are empty """ if points.ndim != 2 or points.shape[1] != 3: raise ValueError(f"Expected points with shape (n, 3), got {points.shape}.") @@ -135,6 +147,11 @@ def calculate_center_point(points: np.ndarray) -> np.ndarray: def find_farthest_point(points: np.ndarray, reference_point: np.ndarray) -> np.ndarray: """ Calculate the farthest point inside a point cloud from a reference point. + + :param points: point coordinates with shape (n, 3) + :param reference_point: point from which distances are calculated + :return: point with the greatest distance from the reference point + :raises ValueError: if the reference point or the param points have an invalid shape or are empty """ if points.ndim != 2 or points.shape[1] != 3: raise ValueError(f"Expected points with shape (n, 3), got {points.shape}.") @@ -156,6 +173,12 @@ def find_farthest_point_vdw( ) -> tuple[np.ndarray, float]: """ Calculate the point, whose fdw-"Bubble" is most distant to a reference point. + + :param points: atom coordinates with shape (n, 3) + :param elements: element symbol mapping to each atom coordinate + :param reference_point: point from which distances are calculated + :return: farthest atom coordinate and its distance including the van der Waals radius + :raises ValueError: if the reference point or the param points have an invalid shape or are empty """ if points.ndim != 2 or points.shape[1] != 3: raise ValueError(f"Expected points with shape (n, 3), got {points.shape}.") @@ -179,6 +202,13 @@ def find_intersecting_spheres( spheres: list[dict], reference_sphere: dict, ) -> list[dict]: + """ + Find spheres that intersect or touch a reference sphere. + + :param spheres: spheres represented by center coordinates and radius + :param reference_sphere: sphere against which intersections are checked + :return: spheres intersecting or touching the reference sphere + """ reference_center = np.array(reference_sphere["center"], dtype=float) reference_radius = float(reference_sphere["radius"]) @@ -201,6 +231,12 @@ def meshes_intersect( ) -> bool: """ Determine whether two meshes intersect or touch. + + :param mesh_a: first mesh + :param mesh_b: second mesh + :param distance_tolerance: maximum distance at which meshes count as touching + (standard is 1e-9 because of floating point error) + :return: whether the meshes intersect or are within the distance tolerance """ from trimesh.collision import CollisionManager @@ -213,6 +249,10 @@ def meshes_intersect( def meshes_distance(mesh_a: trimesh.Trimesh, mesh_b: trimesh.Trimesh) -> float: """ Calculate the minimum euclidean distance between two meshes. + + :param mesh_a: first mesh + :param mesh_b: second mesh + :return: minimum Euclidean distance between the meshes """ from trimesh.collision import CollisionManager From 6b492b16aab85d6abb5a6dbd9f50eb386d47f5c0 Mon Sep 17 00:00:00 2001 From: 3dot141592 Date: Mon, 15 Jun 2026 12:49:13 +0200 Subject: [PATCH 22/22] fix tests --- backend/tests/main/test_views_helper.py | 1 + .../protzilla/utilities/test_ptm_helper.py | 26 +++++++++++++++---- 2 files changed, 22 insertions(+), 5 deletions(-) diff --git a/backend/tests/main/test_views_helper.py b/backend/tests/main/test_views_helper.py index b8e38d9cd..9c03c8cf5 100644 --- a/backend/tests/main/test_views_helper.py +++ b/backend/tests/main/test_views_helper.py @@ -92,6 +92,7 @@ def test_get_all_possible_step_names(): "ArbitraryCSVImport", "CrosslinkingValidationWithAngstromDeviation", "CrosslinkingValidationWithAngstromDeviationForMultimer", + "PtmValidation", } steps = get_all_possible_steps() diff --git a/backend/tests/protzilla/utilities/test_ptm_helper.py b/backend/tests/protzilla/utilities/test_ptm_helper.py index f63e97461..b25685614 100644 --- a/backend/tests/protzilla/utilities/test_ptm_helper.py +++ b/backend/tests/protzilla/utilities/test_ptm_helper.py @@ -14,6 +14,7 @@ def test_get_all_ptm_atoms_with_coordinates(): "_atom_site.label_seq_id": [10, 10, 42, 42], "_atom_site.label_comp_id": ["MET", "MET", "SEP", "SEP"], "_atom_site.label_atom_id": ["N", "CA", "P", "O1P"], + "_atom_site.type_symbol": ["N", "C", "P", "O"], "_atom_site.Cartn_x": [1.0, 2.0, 10.0, 11.0], "_atom_site.Cartn_y": [1.1, 2.1, 10.1, 11.1], "_atom_site.Cartn_z": [1.2, 2.2, 10.2, 11.2], @@ -32,9 +33,21 @@ def test_get_all_ptm_atoms_with_coordinates(): atoms = ptm["atoms"] assert len(atoms) == 2 - assert atoms[0] == {"atom_name": "P", "x": 10.0, "y": 10.1, "z": 10.2} + assert atoms[0] == { + "atom_name": "P", + "element": "P", + "x": 10.0, + "y": 10.1, + "z": 10.2, + } - assert atoms[1] == {"atom_name": "O1P", "x": 11.0, "y": 11.1, "z": 11.2} + assert atoms[1] == { + "atom_name": "O1P", + "element": "O", + "x": 11.0, + "y": 11.1, + "z": 11.2, + } def test_get_center_points_and_radius_for_each_ptm(): @@ -69,6 +82,9 @@ def test_get_center_points_and_radius_for_each_ptm(): ], }, ] + for ptm in ptm_list: + for atom in ptm["atoms"]: + atom["element"] = "C" result_list = get_center_points_and_radius_for_each_ptm(ptm_list) @@ -78,15 +94,15 @@ def test_get_center_points_and_radius_for_each_ptm(): ptm_1 = result_list[0] np.testing.assert_allclose(ptm_1["center_point"], [0.0, 0.0, 0.0], atol=1e-6) - assert pytest.approx(ptm_1["radius"], 1e-6) == 1.0 + assert pytest.approx(ptm_1["radius"], 1e-6) == 2.77 ptm_2 = result_list[1] np.testing.assert_allclose(ptm_2["center_point"], [5.0, 5.0, 5.0], atol=1e-6) - assert pytest.approx(ptm_2["radius"], 1e-6) == 0.0 + assert pytest.approx(ptm_2["radius"], 1e-6) == 1.77 ptm_3 = result_list[2] expected_center = [10.0, 10.0, 10.0] np.testing.assert_allclose(ptm_3["center_point"], expected_center, atol=1e-6) - assert pytest.approx(ptm_3["radius"], 1e-6) == 5.0 + assert pytest.approx(ptm_3["radius"], 1e-6) == 6.77