PCA-Based Registration

Shape-based registration using principal component analysis.

Class Reference

class physiotwin4d.RegisterModelsPCA(pca_template_model, pca_eigenvectors, pca_std_deviations, pca_number_of_modes=0, pca_template_model_point_subsample=4, post_pca_transform=None, fixed_distance_map=None, fixed_model=None, reference_image=None, pca_prior_weight=0.0, symmetric_weight=0.5, log_level=20)[source]

Bases: PhysioTwin4DBase

Register PCA-based shape models to images by minimizing a distance metric.

This class implements a registration pipeline for fitting statistical shape models to patient-specific medical images:

PCA Deformable Registration

  • Optimizes PCA coefficients

  • Model equation: P = template + Σ(b_i * std_i * pca_eigenvector_i)

  • Minimizes the mean distance-to-target at the deformed model points P

Optimization Objective:

fixed_distance_map is zero on the target surface and grows with distance away from it (in mm), so the objective is minimized:

E(b) = (1 - w) * mean_i D(P_i(b))          # model -> target
     +       w * mean_j ||Q_j - P_nn(j)||  # target -> model
     + lambda * Σ b_i²                     # Mahalanobis prior

w is symmetric_weight and lambda is pca_prior_weight. Because the deformation is linear in b, the gradient is analytic and is supplied to the optimizer directly.

Coordinate frames:

The eigenvectors are directions in the statistical model’s own training frame, so they are only valid when added to a template in that same frame. Any rigid/affine alignment to the target must be supplied as post_pca_transform, which is applied after the deformation, rather than pre-applied to pca_template_model.

pca_template_model

Mean shape model

Type:

pv.DataSet

pca_eigenvectors

PCA eigenvectors/components (modes × n_points*3)

Type:

np.ndarray

pca_std_deviations

Standard deviations per mode (modes,)

Type:

np.ndarray

fixed_distance_map

Distance map of the target, in mm

Type:

itk.Image

fixed_model

Target model, when one was supplied. Required for the symmetric (target-to-model) term.

Type:

pv.DataSet

pca_number_of_modes

Number of PCA modes available

Type:

int

registered_model_pca_coefficients

Optimized PCA coefficients

Type:

np.ndarray

registered_model

Final registered and deformed model

Type:

pv.DataSet

post_pca_transform

Transform to apply after PCA registration

Type:

itk.Transform

forward_point_transform

POINT transform mapping template points -> registered/target points; use it to warp the template model/landmarks onto the target. Its orientation is opposite to an image-registration forward_transform (see docs/developer/transform_conventions). Does not include the post-PCA transform.

Type:

itk.DisplacementFieldTransform

inverse_point_transform

POINT transform mapping target points -> template points. Does not include the post-PCA transform.

Type:

itk.DisplacementFieldTransform

Example

>>> # Load PCA model data
>>> pca_template_model = pv.read('pca_All_mean.vtk')
>>> with open('pca.json', 'r') as f:
...     pca_data = json.load(f)
>>> pca_group_data = pca_data['All']
>>> pca_std_deviations = np.sqrt(np.array(pca_group_data['eigenvalues']))
>>> pca_eigenvectors = np.array(pca_group_data['components'])
>>>
>>> # Initialize registrar with loaded data
>>> registrar = RegisterModelsPCA(
...     pca_template_model=pca_template_model,
...     pca_eigenvectors=pca_eigenvectors,
...     pca_std_deviations=pca_std_deviations,
... )
>>>
>>> # Run full registration pipeline
>>> result = registrar.register(pca_number_of_modes=10)
>>>
>>> # Save registered model
>>> result['registered_model'].save('registered_heart.vtk')
>>>
>>> # Print optimization results
>>> print(f'Final mean distance: {result["mean_distance"]:.2f}')
>>> print(f'PCA coefficients: {result["pca_coefficients"]}')
__init__(pca_template_model, pca_eigenvectors, pca_std_deviations, pca_number_of_modes=0, pca_template_model_point_subsample=4, post_pca_transform=None, fixed_distance_map=None, fixed_model=None, reference_image=None, pca_prior_weight=0.0, symmetric_weight=0.5, log_level=20)[source]

Initialize the PCA-based model-to-image registration.

Parameters:

  • pca_template_model (DataSet) – PyVista model containing the mean 3D shape model (unstructured grid or polydata). It must be in the same frame the PCA modes were trained in; supply any alignment to the target as post_pca_transform instead of pre-applying it here.

  • pca_eigenvectors (ndarray) – Numpy array of PCA eigenvectors/components. Shape: (modes, n_points*3) Each row is a flattened eigenmode with 3D displacements: [x1,y1,z1, x2,y2,z2, …]

  • pca_std_deviations (ndarray) – Numpy array of standard deviations per PCA mode. Shape: (modes,) These are the square roots of pca_eigenvalues

  • pca_number_of_modes (int) – Number of PCA modes to use. Default: 0 (use all)

  • pca_template_model_point_subsample (int) – Step size for subsampling model points. Default: 4

  • post_pca_transform (Optional[Transform]) – Optional ITK transform to apply after PCA registration. Default: None

  • fixed_distance_map (Optional[Image]) – ITK image providing the distance map, in mm. Default: None

  • fixed_model (Optional[DataSet]) – PyVista model used to compute the distance map, if one isn’t provided. Also supplies the target points for the symmetric term.

  • reference_image (Optional[Image]) – ITK image providing coordinate frame for computing the distance map.

  • pca_prior_weight (float) – Weight (in mm) of the Mahalanobis shape prior lambda * sum(b_i**2). Because b is expressed in standard deviations, this term is the squared Mahalanobis distance in shape space and makes the fit a MAP estimate rather than a pure data fit constrained only by the coefficient bounds. Default: 0.0 (prior disabled).

  • symmetric_weight (float) – Weight in [0, 1] of the target-to-model distance term. 0.0 measures model-to-target only, which lets the model satisfy the metric while covering just part of the target. Requires fixed_model; ignored with a warning when only a distance map is available. Default: 0.5

  • log_level (int | str) – Logging level (logging.DEBUG, logging.INFO, logging.WARNING). Default: logging.INFO

Raises:

ValueError – If pca_eigenvector dimensions don’t match model points, if the mode counts disagree, or if neither a distance map nor a fixed model plus reference image is provided.

classmethod from_json(pca_template_model, pca_json_filename, pca_number_of_modes=0, pca_template_model_point_subsample=4, post_pca_transform=None, fixed_distance_map=None, fixed_model=None, reference_image=None, pca_prior_weight=0.0, symmetric_weight=0.5, log_level=20)[source]

Create RegisterModelsPCA from PCA model JSON file.

This method reads PCA statistical shape model data from a JSON file containing eigenvalues and principal component vectors.

The JSON file must contain: - ‘eigenvalues’: Array of eigenvalues (variance) for each component - ‘components’: Array of principal component vectors (flattened shape deformations)

Parameters:

  • pca_template_model (DataSet) – Mean surface mesh to use as template

  • pca_json_filename (str) – Path to the PCA model JSON file

  • pca_number_of_modes (int) – Number of PCA modes to use. Default: 0 (use all)

  • pca_template_model_point_subsample (int) – Step size for subsampling model points. Default: 4

  • post_pca_transform (Optional[Transform]) – Optional ITK transform to apply after PCA registration. Default: None

  • fixed_distance_map (Optional[Image]) – ITK image providing the distance values for registration. If None, must be set later before registration.

  • fixed_model (Optional[DataSet]) – Target surface mesh to register to. Default: None

  • reference_image (Optional[Image]) – Reference image defining coordinate space. Default: None

  • pca_prior_weight (float) – Weight (mm) of the Mahalanobis shape prior. Default: 0.0

  • symmetric_weight (float) – Weight of the target-to-model term. Default: 0.5

  • log_level (int | str) – Logging level (logging.DEBUG, logging.INFO, logging.WARNING). Default: logging.INFO

Return type:

Self

Returns:

RegisterModelsPCA instance

Raises:

Example

>>> registrar = RegisterModelsPCA.from_json(
...     pca_template_model=pca_template_model,
...     pca_json_filename='path/to/pca_model.json',
...     fixed_model=fixed_model,
...     reference_image=reference_image,
... )
classmethod from_pca_model(pca_template_model, pca_model, pca_number_of_modes=0, pca_template_model_point_subsample=4, post_pca_transform=None, fixed_distance_map=None, fixed_model=None, reference_image=None, pca_prior_weight=0.0, symmetric_weight=0.5, log_level=20)[source]

Create RegisterModelsPCA from a PCA model dictionary.

The dict must match the structure produced by WorkflowCreateStatisticalModel (key pca_model): explained_variance_ratio, eigenvalues, components.

Parameters:

  • pca_template_model (DataSet) – Mean surface mesh to use as template

  • pca_model (dict) – PCA model dict with ‘eigenvalues’ and ‘components’ (and optionally ‘explained_variance_ratio’)

  • pca_number_of_modes (int) – Number of PCA modes to use. Default: 0 (use all)

  • pca_template_model_point_subsample (int) – Step size for subsampling model points. Default: 4

  • post_pca_transform (Optional[Transform]) – Optional ITK transform to apply after PCA registration.

  • fixed_distance_map (Optional[Image]) – ITK image providing the distance values for registration.

  • fixed_model (Optional[DataSet]) – Target surface mesh to register to.

  • reference_image (Optional[Image]) – Reference image defining coordinate space.

  • pca_prior_weight (float) – Weight (mm) of the Mahalanobis shape prior. Default: 0.0

  • symmetric_weight (float) – Weight of the target-to-model term. Default: 0.5

  • log_level (int | str) – Logging level.

Return type:

Self

Returns:

RegisterModelsPCA instance

Raises:

ValueError – If required keys are missing or dimensions invalid

set_fixed_model(fixed_model, reference_image)[source]

Set the fixed model for registration and rebuild its distance map.

Parameters:

  • fixed_model (UnstructuredGrid) – PyVista model used to compute the distance map.

  • reference_image (Optional[Image]) – ITK image providing coordinate frame for computing the distance map.

Return type:

None

set_fixed_distance_map(fixed_distance_map)[source]

Set the distance map used as the registration target.

Parameters:

fixed_distance_map (Optional[Image]) – ITK image providing distance data, in mm

Return type:

None

set_pca_template_model(pca_template_model)[source]

Set the average model for registration.

Parameters:

pca_template_model (UnstructuredGrid) – PyVista model containing the mean 3D shape model (unstructured grid or polydata)

Return type:

None

transform_template_model()[source]

Create the final registered model by applying PCA deformation.

Return type:

DataSet

Returns:

Final registered and deformed model as a PyVista dataset.

Raises:

ValueError – If registration has not been performed

transform_point(point, include_post_pca_transform=True)[source]

Transform an arbitrary point through the PCA deformation field.

Parameters:

  • point (Point) – ITK point to transform (itk.Point[itk.D, 3])

  • include_post_pca_transform (bool) – Also apply post_pca_transform. Default: True

Return type:

Point

Returns:

Transformed ITK point

Raises:

ValueError – If compute_pca_transforms() has not been called yet

Notes

This samples the approximated deformation field built by compute_pca_transforms(), which is splatted and blurred, so it does not reproduce transform_template_model() exactly; the RMS of that difference is logged when the field is built. Points outside the field’s reference image are not displaced.

Example

>>> p = itk.Point[itk.D, 3]()
>>> p[0], p[1], p[2] = 10.0, 20.0, 30.0
>>> transformed_p = registrar.transform_point(p)
compute_pca_transforms(reference_image, blur_sigma=2.5)[source]

Compute PCA transforms.

The field is built by splatting the per-point PCA displacements onto the reference grid and blurring them, so it only approximates the exact per-point deformation. The RMS of that approximation error, and of the forward/inverse round trip, are both logged.

Parameters:

  • reference_image (Image) – ITK image providing the coordinate frame for the field.

  • blur_sigma (float) – Sigma for Gaussian blurring of the deformation field. Default: 2.5

Returns:

  • ‘forward_point_transform’: POINT transform mapping template points -> target points (warps the template onto the target)

  • ’inverse_point_transform’: POINT transform mapping target points -> template points

Return type:

Dictionary containing

Note

These are point transforms, oriented opposite to image-registration transforms; see docs/developer/transform_conventions. Neither includes post_pca_transform.

classmethod get_log_classes()

Get the list of classes currently showing logs.

Return type:

list[str]

Returns:

List of class names that are allowed to show logs. Empty list if filter is disabled (all classes shown).

Example

>>> classes = PhysioTwin4DBase.get_log_classes()
>>> print(classes)
['RegisterModelsPCA', 'WorkflowFitStatisticalModelToPatient']
log_critical(message, *args)

Log a critical message with optional %-style formatting.

Parameters:

  • message (str) – The critical message to log (can contain %-style placeholders)

  • *args (Any) – Arguments for %-style string formatting

Return type:

None

Example

>>> self.log_critical('System failure at %s', timestamp)
>>> self.log_critical('Critical error: %(msg)s', {'msg': 'Out of memory'})
log_debug(message, *args)

Log a debug message with optional %-style formatting.

Parameters:

  • message (str) – The debug message to log (can contain %-style placeholders)

  • *args (Any) – Arguments for %-style string formatting

Return type:

None

Example

>>> self.log_debug('Processing %s with %d items', filename, count)
>>> self.log_debug('Value is %(value)d', {'value': 42})
log_error(message, *args)

Log an error message with optional %-style formatting.

Parameters:

  • message (str) – The error message to log (can contain %-style placeholders)

  • *args (Any) – Arguments for %-style string formatting

Return type:

None

Example

>>> self.log_error('Failed to load %s: %s', filename, error_msg)
>>> self.log_error('Error code: %(code)d', {'code': 404})
log_info(message, *args)

Log an info message with optional %-style formatting.

Parameters:

  • message (str) – The info message to log (can contain %-style placeholders)

  • *args (Any) – Arguments for %-style string formatting

Return type:

None

Example

>>> self.log_info('Loading file: %s', filepath)
>>> self.log_info('Iteration %(iter)d of %(total)d', {'iter': 5, 'total': 10})
log_progress(current, total, prefix='Progress')

Log progress information.

Parameters:

  • current (int) – Current step/iteration number

  • total (int) – Total number of steps/iterations

  • prefix (str) – Prefix text for the progress message. Default: ‘Progress’

Return type:

None

Example

>>> for i in range(100):
...     self.log_progress(i + 1, 100)
>>> self.log_progress(5, 10, prefix='Processing')

Note

For custom formatted progress messages, use log_info() directly: >>> self.log_info(‘Loading %s: %d/%d’, filename, current, total)

log_section(title, *args, width=70, char='=')

Log a formatted section header with optional %-style formatting.

Useful for visually separating major sections of output.

Parameters:

  • title (str) – The section title (can contain %-style placeholders)

  • *args (Any) – Arguments for %-style string formatting of title

  • width (int) – Total width of the header line. Default: 70

  • char (str) – Character to use for the header line. Default: ‘=’

Return type:

None

Example

>>> self.log_section('Stage 1: Initialization')
>>> self.log_section('Processing file: %s', filename)
>>> self.log_section('Stage %(num)d: %(name)s', {'num': 2, 'name': 'Analysis'})
# Outputs:
# ======================================================================
# Stage 2: Analysis
# ======================================================================
log_warning(message, *args)

Log a warning message with optional %-style formatting.

Parameters:

  • message (str) – The warning message to log (can contain %-style placeholders)

  • *args (Any) – Arguments for %-style string formatting

Return type:

None

Example

>>> self.log_warning('Memory usage at %d%%', usage_percent)
>>> self.log_warning('Parameter %(name)s out of range', {'name': 'threshold'})
classmethod set_log_all_classes()

Enable logging output from all PhysioTwin4D classes.

Disables the class filter so all classes show their logs.

Example

>>> PhysioTwin4DBase.set_log_all_classes()
>>> # Now all classes will show their logs
Return type:

None

classmethod set_log_classes(class_names)

Set which classes should show their logging output.

Only log messages from the specified classes will be displayed. All other classes will have their logs hidden.

Parameters:

class_names (list[str]) – List of class names to show logs from. Example: [“RegisterModelsPCA”, “WorkflowFitStatisticalModelToPatient”]

Return type:

None

Example

>>> PhysioTwin4DBase.set_log_classes(['RegisterModelsPCA'])
>>> # Now only RegisterModelsPCA logs will be shown
classmethod set_log_level(log_level)

Set the logging level for all PhysioTwin4D classes.

Parameters:

log_level (int | str) – Logging level. Can be an integer (logging.DEBUG, logging.INFO, logging.WARNING, logging.ERROR, logging.CRITICAL) or a string (‘DEBUG’, ‘INFO’, ‘WARNING’, ‘ERROR’, ‘CRITICAL’).

Return type:

None

Example

>>> import logging
>>> PhysioTwin4DBase.set_log_level(logging.DEBUG)
>>> # or
>>> PhysioTwin4DBase.set_log_level('DEBUG')
register(pca_number_of_modes=0, pca_coefficient_bounds=3.5, method='L-BFGS-B', max_iterations=100)[source]

Optimize PCA coefficients to deform the model onto the target.

Parameters:

  • pca_number_of_modes (int) – Number of PCA modes to use. Default: 0 (use all available modes)

  • pca_coefficient_bounds (float) – PCA coefficient bounds (±std devs). Default: 3.5

  • method (str) – Optimization method for scipy.optimize.minimize. Default: ‘L-BFGS-B’ (supports bounds)

  • max_iterations (int) – Maximum number of optimization iterations. Default: 100

Returns:

  • ‘registered_model’: Final registered PyVista model

  • ’pca_coefficients’: Optimized PCA coefficients

  • ’mean_distance’: Final objective value, in mm

Return type:

Dictionary containing

Raises:

ValueError – If the distance map is not set

Example

>>> result = registrar.register(pca_number_of_modes=10)
>>> result['registered_model'].save('registered_heart.vtk')

Navigation

Distance Map Registration | Model Registration Modules | USD Generation Modules