diff --git a/docs/monai_lib.rst b/docs/monai_lib.rst
new file mode 100644
index 00000000..32b5dd43
--- /dev/null
+++ b/docs/monai_lib.rst
@@ -0,0 +1,50 @@
+.. _monai_lib:
+
+MONAI
+=====
+
+`MONAI`_ is a widely-used library for training and evaluating models on medical imaging
+data. In order to use MONAI with highdicom, the ``monai`` python package must be
+installed separately. Version 1.3.0 or later is required.
+
+.. _monai_vol:
+
+Volume Conversions
+------------------
+
+Highdicom supports conversions with the ``monai.data.MetaTensor`` class through the
+:meth:`highdicom.Volume.to_monai` and :meth:`highdicom.Volume.from_monai` methods.
+Unlike highdicom, which exclusively uses the "LPS" convention, MONAI supports both
+"RAS" and "LPS" MetaTensors. By default, highdicom converts to the "RAS" convention.
+
+Creating an MONAI MetaTensor from a Volume:
+
+.. code-block:: python
+
+    import highdicom as hd
+
+
+    vol = hd.Volume(...)
+
+    metatensor = vol.to_monai()
+
+Creating a volume from a MONAI MetaTensor:
+
+.. code-block:: python
+
+    import monai
+    import highdicom as hd
+
+
+    metatensor = monai.data.MetaTensor(...)
+
+    vol = hd.Volume.from_monai(
+        metatensor=metatensor,
+        coordinate_system='PATIENT',
+        frame_of_reference_uid=None,
+        channels=None,
+        channel_dim=None
+    )
+
+
+.. _`MONAI`: https://project-monai.github.io/
diff --git a/docs/other_libraries.rst b/docs/other_libraries.rst
index 566dabd8..5a6ddb70 100644
--- a/docs/other_libraries.rst
+++ b/docs/other_libraries.rst
@@ -23,3 +23,4 @@ library enables.
    itk_lib
    sitk_lib
    nib_lib
+   monai_lib
diff --git a/pyproject.toml b/pyproject.toml
index 9e45200f..7de06911 100644
--- a/pyproject.toml
+++ b/pyproject.toml
@@ -51,6 +51,7 @@ libjpeg = [
 itk = ["itk>=5.4.0"]
 simpleitk = ["SimpleITK>=2.2.1"]
 nibabel = ["nibabel>=4.0.0"]
+monai = ["monai>=1.3.0"]
 
 [dependency-groups]
 test = [
diff --git a/src/highdicom/volume.py b/src/highdicom/volume.py
index ac8f294f..fb9c7bd6 100644
--- a/src/highdicom/volume.py
+++ b/src/highdicom/volume.py
@@ -4241,6 +4241,228 @@ def from_nibabel(
             from_reference_convention='RAS'
         )
 
+    def to_monai(
+        self,
+        convert_to_ras: bool = True,
+        squeeze_dim: None | int = None,
+        ensure_channel_first: bool = False,
+    ) -> 'monai.data.MetaTensor':  # noqa: F821
+        """Convert the volume to a ``monai.data.MetaTensor``.
+
+        This method requires an optional dependency to be installed
+        separately from highdicom, specifically ``monai``.
+
+        The Volume is converted to a ``monai.data.Metatensor``. Spatial metadata
+        is preserved through the affine array. By default, metatensors use "RAS"
+        convention, however they can support "LPS" similar to highdicom if
+        specified.
+
+        Parameters
+        ----------
+        convert_to_ras: bool
+            Whether to convert the affine matrix from 'LPS' to 'RAS' convention.
+        squeeze_dim: None | int
+            Index of a singleton dimension to squeeze (must be a spatial
+            dimension). Defaults to None.
+        ensure_channel_first: bool
+            Whether to convert to a channel first metatensor. Defaults to False.
+
+        Returns
+        -------
+        monai.data.MetaTensor:
+            MetaTensor constructed from the volume.
+
+        Raises
+        ------
+        ValueError
+            When the volume has more than one channel dimension.
+        ValueError
+            When squeeze_dim does not correspond to a spatial dimension.
+
+        """
+        func = self.to_monai
+        monai = import_optional_dependency(
+            module_name='monai',
+            feature=f'{func.__module__}.{func.__qualname__}'
+        )
+
+        ImageStatsKeys = monai.utils.enums.ImageStatsKeys
+        MetaKeys = monai.utils.enums.MetaKeys
+
+        if self.number_of_channel_dimensions > 1:
+            raise ValueError(
+                'Monai conversion does not currently support'
+                ' volumes with multiple channel dimensions.'
+            )
+
+        meta = {}
+        if convert_to_ras:
+            space = monai.utils.enums.SpaceKeys.RAS
+
+        else:
+            space = monai.utils.enums.SpaceKeys.LPS
+
+        affine = self.get_affine(space.value)
+        array = np.ascontiguousarray(self.array)
+
+        keep_cols = [i for i in range(3) if i != squeeze_dim]
+        if squeeze_dim is not None:
+            if squeeze_dim not in [0, 1, 2]:
+                raise ValueError(
+                    'If provided, `squeeze_dim` must be a spatial'
+                    ' dimension (0, 1, 2).'
+                )
+
+            if array.shape[squeeze_dim] != 1:
+                raise ValueError(
+                    f'`squeeze_dim={squeeze_dim}` does not correspond'
+                    ' to a singleton dimension. Array has shape:'
+                    f' {array.shape}.'
+                )
+
+            u = affine[:3, keep_cols[0]]
+            v = affine[:3, keep_cols[1]]
+            origin = affine[:3, 3]
+
+            # define new orthonormal basis (e_0, e_1)
+            e_0 = u / np.linalg.norm(u)
+            v_perp = v - (np.dot(v, e_0) * e_0)
+            e_1 = v_perp / np.linalg.norm(v_perp)
+
+            affine = np.array(
+                [[np.dot(u, e_0), np.dot(v, e_0), np.dot(origin, e_0)],
+                 [np.dot(u, e_1), np.dot(v, e_1), np.dot(origin, e_1)],
+                 [0, 0, 1]]
+            )
+            array = array.squeeze(squeeze_dim)
+
+        meta[MetaKeys.SPACE] = space
+        meta[ImageStatsKeys.SPACING] = tuple([
+            self.spacing[i] for i in keep_cols
+        ])
+        meta[MetaKeys.SPATIAL_SHAPE] = np.array(self.spatial_shape)[keep_cols]
+        meta[MetaKeys.ORIGINAL_AFFINE] = affine.copy()
+        meta[MetaKeys.AFFINE] = affine.copy()
+        meta[MetaKeys.ORIGINAL_CHANNEL_DIM] = (
+            -1 if self.ndim > 3 else float("nan")
+        )
+
+        metatensor = monai.data.MetaTensor(array, meta=meta)
+
+        if ensure_channel_first:
+            metatensor = monai.transforms.EnsureChannelFirst()(metatensor)
+
+        return metatensor
+
+    @classmethod
+    def from_monai(
+        cls,
+        metatensor: 'monai.data.MetaTensor',  # noqa: F821
+        coordinate_system: CoordinateSystemNames | str = 'PATIENT',
+        frame_of_reference_uid: str | None = None,
+        channels: dict[
+            BaseTag | int | str | ChannelDescriptor,
+            Sequence[int | str | float | Enum]
+        ] | None = None,
+        channel_dim: int = 0
+    ) -> Self:
+        """Construct a Volume from a ``monai.data.MetaTensor``.
+
+        This method requires an optional dependency to be installed
+        separately from highdicom, specifically ``monai``.
+
+        The ``monai.data.MetaTensor`` is converted to a 3D Volume.
+        Spatial metadata is preserved through the affine array. By default,
+        metatensors use "RAS" convention, however they can support "LPS" similar
+        to highdicom. For correct conversions, ensure that the metatensor's
+        metadata specifies the appropriate convention through the space key
+        (``monai.utils.enums.SPACE``).
+
+        Parameters
+        ----------
+        metatensor: monai.data.MetaTensor
+            A ``monai.data.MetaTensor`` to convert to a volume.
+        coordinate_system: highdicom.CoordinateSystemNames | str
+            Coordinate system (``"PATIENT"`` or ``"SLIDE"``) in which the volume
+            is defined.
+        frame_of_reference_uid: Union[str, None], optional
+            Frame of reference UID for the frame of reference, if known.
+        channels: dict[int | str | ChannelDescriptor, Sequence[int | str | float | Enum]] | None, optional
+            Specification of channels of the array. Channels are additional
+            dimensions of the array beyond the three spatial dimensions. For
+            each such additional dimension (if any), an item in this dictionary
+            is required to specify the meaning. The dictionary key specifies
+            the meaning of the dimension, which must be either an instance of
+            highdicom.ChannelDescriptor, specifying a DICOM tag whose attribute
+            describes the channel, a a DICOM keyword describing a DICOM
+            attribute, or an integer representing the tag of a DICOM attribute.
+            The corresponding item of the dictionary is a sequence giving the
+            value of the relevant attribute at each index in the array. The
+            insertion order of the dictionary is significant as it is used to
+            match items to the corresponding dimensions of the array (the first
+            item in the dictionary corresponds to axis 3 of the array and so
+            on).
+        channel_dim: int
+            Channel dimension index of the metatensor. Defaults to 0.
+
+        Returns
+        -------
+        highdicom.Volume:
+            Volume constructed from the monai image.
+
+        Raises
+        ------
+        ValueError
+            When the metatensor has more than 4 dimensions (multiple channel
+            dimensions are unsupported).
+        ValueError
+            When there are multiple channels and the `channels` argument is
+            not provided.
+
+        """  # noqa: E501
+        func = cls.from_monai
+        monai = import_optional_dependency(
+            module_name='monai',
+            feature=f'{func.__module__}.{func.__qualname__}'
+        )
+
+        if metatensor.ndim > 3:
+            if metatensor.ndim > 4:
+                raise ValueError(
+                    'Monai conversion does not currently support'
+                    ' volumes with multiple channel dimensions.'
+                )
+
+            channel_dim = channel_dim % metatensor.ndim
+            channel_last_perm = [
+                i for i in range(metatensor.ndim) if i != channel_dim
+            ] + [channel_dim]
+            metatensor = metatensor.permute(channel_last_perm)
+
+            if channels is None:
+                if metatensor.shape[-1] == 1:
+                    metatensor = metatensor.squeeze(-1)
+
+                else:
+                    raise ValueError(
+                        'Monai conversion requires `channels` be specified'
+                        ' for volumes with >=2 channels.'
+                    )
+
+        array = metatensor.cpu().detach().numpy()
+
+        return cls(
+            array=array,
+            affine=metatensor.affine.cpu().detach().numpy(),
+            coordinate_system=coordinate_system,
+            frame_of_reference_uid=frame_of_reference_uid,
+            channels=channels,
+            from_reference_convention=metatensor.meta.get(
+                monai.utils.enums.MetaKeys.SPACE,
+                monai.utils.enums.SpaceKeys.RAS
+            ).value
+        )
+
 
 class VolumeToVolumeTransformer:
 
diff --git a/tests/test_monai.py b/tests/test_monai.py
new file mode 100644
index 00000000..b385bfab
--- /dev/null
+++ b/tests/test_monai.py
@@ -0,0 +1,975 @@
+import numpy as np
+import tempfile
+import pytest
+import pydicom
+import zipfile
+import re
+
+from pathlib import Path
+from typing import Sequence
+from highdicom import Volume, ChannelDescriptor, get_volume_from_series
+from highdicom.spatial import (
+    get_closest_patient_orientation,
+    convert_affine_to_convention
+)
+from highdicom._dependency_utils import import_optional_dependency
+from highdicom.seg import segread
+from .utils import (
+    DCM_QA_MPRAGE,
+    DCM_QA_ME,
+    DCM_QA_PDT2,
+    TEST_DATA,
+    read_multiframe_ct_volume,
+    read_ct_series_volume,
+    urldownload_with_retry
+)
+
+
+try:
+    monai = import_optional_dependency('monai', feature='monai tests')
+    ImageStatsKeys = monai.utils.enums.ImageStatsKeys
+    MetaKeys = monai.utils.enums.MetaKeys
+
+except Exception:
+    pytest.skip("Optional dependency not available", allow_module_level=True)
+
+
+def read_github_zip_volume_and_metatensor(url: str):
+    with tempfile.TemporaryDirectory() as temp_dir:
+        zipfilename = Path(temp_dir) / Path(url).name
+        urldownload_with_retry(url, zipfilename)
+
+        with zipfile.ZipFile(zipfilename, 'r') as zf:
+            zf.extractall(temp_dir)
+
+        series = [pydicom.dcmread(f) for f in Path(temp_dir).glob('**/*.dcm')]
+
+        try:
+            pydicom.config.enforce_valid_values = False
+            metatensor = monai.transforms.LoadImage(reader="PydicomReader")(
+                Path(str(zipfilename)[:-4])
+            )
+
+        finally:
+            pydicom.config.enforce_valid_values = True
+
+    return get_volume_from_series(series), series, metatensor
+
+
+def read_github_series_volume_and_metatensor(urls: Sequence[str]):
+    series = []
+    with tempfile.TemporaryDirectory() as temp_dir:
+        for url in urls:
+            filename = Path(temp_dir) / Path(url).name
+            urldownload_with_retry(url, filename)
+
+            series.append(pydicom.dcmread(filename))
+
+        try:
+            pydicom.config.enforce_valid_values = False
+            metatensor = monai.transforms.LoadImage(
+                reader="PydicomReader"
+            )(temp_dir)
+
+        finally:
+            pydicom.config.enforce_valid_values = True
+
+    return get_volume_from_series(series), series, metatensor
+
+
+def validate_affine(
+    affine_3d: np.ndarray,
+    affine_2d: np.ndarray,
+    squeeze_dim: int,
+    shape: Sequence[int]
+):
+    keep_cols = [i for i in range(3) if i != squeeze_dim]
+    u = affine_3d[:3, keep_cols[0]]
+    v = affine_3d[:3, keep_cols[1]]
+    origin = affine_3d[:3, 3]
+
+    e_0 = u / np.linalg.norm(u)
+    v_perp = v - (np.dot(v, e_0) * e_0)
+    e_1 = v_perp / np.linalg.norm(v_perp)
+
+    basis = np.column_stack([e_0, e_1])
+
+    origin_2d = basis.T @ origin
+    origin_offset = origin - basis @ origin_2d
+
+    for x in range(shape[keep_cols[0]]):
+        for y in range(shape[keep_cols[1]]):
+            vox_3d = np.zeros(3)
+            vox_3d[keep_cols[0]] = x
+            vox_3d[keep_cols[1]] = y
+            world_3d = affine_3d[:3, :3] @ vox_3d + affine_3d[:3, 3]
+
+            vox_2d = np.array([x, y])
+            world_2d = affine_2d[:2, :2] @ vox_2d + affine_2d[:2, 2]
+            world_3d_recon = basis @ world_2d + origin_offset
+
+            assert np.allclose(
+                world_3d,
+                world_3d_recon
+            )
+
+
+@pytest.mark.parametrize(
+    'vol',
+    [
+        # testdata_files
+        read_multiframe_ct_volume(),
+        read_ct_series_volume(),
+        # different orientations
+        Volume.from_attributes(
+            array=np.random.randint(0, 100, (16, 16, 16)),
+            image_position=[0.0, 0.0, 0.0],
+            image_orientation=[1.0, 0.0, 0.0, 0.0, 1.0, 0.0],
+            pixel_spacing=[1.0, 1.0],
+            spacing_between_slices=1.0,
+            coordinate_system='PATIENT'
+        ).to_patient_orientation('RAF'),
+        Volume.from_attributes(
+            array=np.random.randint(0, 100, (16, 16, 16)),
+            image_position=[0.0, 0.0, 0.0],
+            image_orientation=[1.0, 0.0, 0.0, 0.0, 1.0, 0.0],
+            pixel_spacing=[1.0, 1.0],
+            spacing_between_slices=1.0,
+            coordinate_system='PATIENT'
+        ).to_patient_orientation('RAH'),
+        Volume.from_attributes(
+            array=np.random.randint(0, 100, (16, 16, 16)),
+            image_position=[0.0, 0.0, 0.0],
+            image_orientation=[1.0, 0.0, 0.0, 0.0, 1.0, 0.0],
+            pixel_spacing=[1.0, 1.0],
+            spacing_between_slices=1.0,
+            coordinate_system='PATIENT'
+        ).to_patient_orientation('RPF'),
+        Volume.from_attributes(
+            array=np.random.randint(0, 100, (16, 16, 16)),
+            image_position=[0.0, 0.0, 0.0],
+            image_orientation=[1.0, 0.0, 0.0, 0.0, 1.0, 0.0],
+            pixel_spacing=[1.0, 1.0],
+            spacing_between_slices=1.0,
+            coordinate_system='PATIENT'
+        ).to_patient_orientation('RPH'),
+        Volume.from_attributes(
+            array=np.random.randint(0, 100, (16, 16, 16)),
+            image_position=[0.0, 0.0, 0.0],
+            image_orientation=[1.0, 0.0, 0.0, 0.0, 1.0, 0.0],
+            pixel_spacing=[1.0, 1.0],
+            spacing_between_slices=1.0,
+            coordinate_system='PATIENT'
+        ).to_patient_orientation('LAF'),
+        Volume.from_attributes(
+            array=np.random.randint(0, 100, (16, 16, 16)),
+            image_position=[0.0, 0.0, 0.0],
+            image_orientation=[1.0, 0.0, 0.0, 0.0, 1.0, 0.0],
+            pixel_spacing=[1.0, 1.0],
+            spacing_between_slices=1.0,
+            coordinate_system='PATIENT'
+        ).to_patient_orientation('LAH'),
+        Volume.from_attributes(
+            array=np.random.randint(0, 100, (16, 16, 16)),
+            image_position=[0.0, 0.0, 0.0],
+            image_orientation=[1.0, 0.0, 0.0, 0.0, 1.0, 0.0],
+            pixel_spacing=[1.0, 1.0],
+            spacing_between_slices=1.0,
+            coordinate_system='PATIENT'
+        ).to_patient_orientation('LPF'),
+        Volume.from_attributes(
+            array=np.random.randint(0, 100, (16, 16, 16)),
+            image_position=[0.0, 0.0, 0.0],
+            image_orientation=[1.0, 0.0, 0.0, 0.0, 1.0, 0.0],
+            pixel_spacing=[1.0, 1.0],
+            spacing_between_slices=1.0,
+            coordinate_system='PATIENT'
+        ).to_patient_orientation('LPH'),
+        Volume.from_attributes(
+            array=np.random.randint(0, 100, (16, 16, 16)),
+            image_position=[0.0, 0.0, 0.0],
+            image_orientation=[1.0, 0.0, 0.0, 0.0, 1.0, 0.0],
+            pixel_spacing=[1.0, 1.0],
+            spacing_between_slices=1.0,
+            coordinate_system='PATIENT'
+        ).to_patient_orientation('HLP'),
+        Volume.from_attributes(
+            array=np.random.randint(0, 100, (16, 16, 16)),
+            image_position=[0.0, 0.0, 0.0],
+            image_orientation=[1.0, 0.0, 0.0, 0.0, 1.0, 0.0],
+            pixel_spacing=[1.0, 1.0],
+            spacing_between_slices=1.0,
+            coordinate_system='PATIENT'
+        ).to_patient_orientation('FPR'),
+        Volume.from_attributes(
+            array=np.random.randint(0, 100, (16, 16, 16)),
+            image_position=[0.0, 0.0, 0.0],
+            image_orientation=[1.0, 0.0, 0.0, 0.0, 1.0, 0.0],
+            pixel_spacing=[1.0, 1.0],
+            spacing_between_slices=1.0,
+            coordinate_system='PATIENT'
+        ).to_patient_orientation('HRP'),
+        # isotropic
+        Volume.from_attributes(
+            array=np.random.randint(0, 100, (16, 16, 16)),
+            image_position=[0.0, 0.0, 0.0],
+            image_orientation=[1.0, 0.0, 0.0, 0.0, 1.0, 0.0],
+            pixel_spacing=[0.5, 0.5],
+            spacing_between_slices=0.5,
+            coordinate_system='PATIENT'
+        ),
+        Volume.from_attributes(
+            array=np.random.randint(0, 100, (16, 16, 16)),
+            image_position=[0.0, 0.0, 0.0],
+            image_orientation=[1.0, 0.0, 0.0, 0.0, 1.0, 0.0],
+            pixel_spacing=[2.0, 2.0],
+            spacing_between_slices=2.0,
+            coordinate_system='PATIENT'
+        ),
+        # anisotropic
+        Volume.from_attributes(
+            array=np.random.randint(0, 100, (16, 16, 16)),
+            image_position=[0.0, 0.0, 0.0],
+            image_orientation=[1.0, 0.0, 0.0, 0.0, 1.0, 0.0],
+            pixel_spacing=[0.5, 0.5],
+            spacing_between_slices=2.0,
+            coordinate_system='PATIENT'
+        ),
+        Volume.from_attributes(
+            array=np.random.randint(0, 100, (16, 16, 16)),
+            image_position=[0.0, 0.0, 0.0],
+            image_orientation=[1.0, 0.0, 0.0, 0.0, 1.0, 0.0],
+            pixel_spacing=[2.0, 0.5],
+            spacing_between_slices=0.5,
+            coordinate_system='PATIENT'
+        ),
+        Volume.from_attributes(
+            array=np.random.randint(0, 100, (16, 16, 16)),
+            image_position=[0.0, 0.0, 0.0],
+            image_orientation=[1.0, 0.0, 0.0, 0.0, 1.0, 0.0],
+            pixel_spacing=[0.5, 2.0],
+            spacing_between_slices=0.5,
+            coordinate_system='PATIENT'
+        ),
+        # non-square
+        Volume.from_attributes(
+            array=np.random.randint(0, 100, (16, 32, 16)),
+            image_position=[0.0, 0.0, 0.0],
+            image_orientation=[1.0, 0.0, 0.0, 0.0, 1.0, 0.0],
+            pixel_spacing=[1.0, 1.0],
+            spacing_between_slices=1.0,
+            coordinate_system='PATIENT'
+        ),
+        Volume.from_attributes(
+            array=np.random.randint(0, 100, (64, 128, 16)),
+            image_position=[0.0, 0.0, 0.0],
+            image_orientation=[1.0, 0.0, 0.0, 0.0, 1.0, 0.0],
+            pixel_spacing=[1.0, 1.0],
+            spacing_between_slices=1.0,
+            coordinate_system='PATIENT'
+        ),
+        # single-slice
+        Volume.from_attributes(
+            array=np.random.randint(0, 100, (16, 16, 1)),
+            image_position=[0.0, 0.0, 0.0],
+            image_orientation=[1.0, 0.0, 0.0, 0.0, 1.0, 0.0],
+            pixel_spacing=[1.0, 1.0],
+            spacing_between_slices=1.0,
+            coordinate_system='PATIENT'
+        ),
+        Volume.from_attributes(
+            array=np.random.randint(0, 100, (16, 1, 16)),
+            image_position=[0.0, 0.0, 0.0],
+            image_orientation=[1.0, 0.0, 0.0, 0.0, 1.0, 0.0],
+            pixel_spacing=[1.0, 1.0],
+            spacing_between_slices=1.0,
+            coordinate_system='PATIENT'
+        ),
+        Volume.from_attributes(
+            array=np.random.randint(0, 100, (1, 16, 16)),
+            image_position=[0.0, 0.0, 0.0],
+            image_orientation=[1.0, 0.0, 0.0, 0.0, 1.0, 0.0],
+            pixel_spacing=[1.0, 1.0],
+            spacing_between_slices=1.0,
+            coordinate_system='PATIENT'
+        ),
+        # random position offset
+        Volume.from_attributes(
+            array=np.random.randint(0, 100, (16, 16, 16)),
+            image_position=[84.40363858, 105.04467386, 143.73326388],
+            image_orientation=[1.0, 0.0, 0.0, 0.0, 1.0, 0.0],
+            pixel_spacing=[1.0, 1.0],
+            spacing_between_slices=1.0,
+            coordinate_system='PATIENT'
+        ),
+        Volume.from_attributes(
+            array=np.random.randint(0, 100, (16, 16, 16)),
+            image_position=[-21.03512292, 35.19549233, -184.42393696],
+            image_orientation=[1.0, 0.0, 0.0, 0.0, 1.0, 0.0],
+            pixel_spacing=[1.0, 1.0],
+            spacing_between_slices=1.0,
+            coordinate_system='PATIENT'
+        ),
+        Volume.from_attributes(
+            array=np.random.randint(0, 100, (16, 16, 16)),
+            image_position=[-197.36060235, 86.22231644, -14.79874245],
+            image_orientation=[1.0, 0.0, 0.0, 0.0, 1.0, 0.0],
+            pixel_spacing=[1.0, 1.0],
+            spacing_between_slices=1.0,
+            coordinate_system='PATIENT'
+        ),
+        # random orientation
+        Volume.from_attributes(
+            array=np.random.randint(0, 100, (16, 16, 16)),
+            image_position=[0.0, 0.0, 0.0],
+            image_orientation=[
+                -0.9267662161157189,
+                0.32283606313442387,
+                -0.1920449348627007,
+                -0.3751482085550474,
+                -0.7693372329674889,
+                0.5170919101937937
+            ],
+            pixel_spacing=[1.0, 1.0],
+            spacing_between_slices=1.0,
+            coordinate_system='PATIENT'
+        ),
+        Volume.from_attributes(
+            array=np.random.randint(0, 100, (16, 16, 16)),
+            image_position=[0.0, 0.0, 0.0],
+            image_orientation=[
+                -0.16411694786392106,
+                -0.7887835415036736,
+                0.5923564400567902,
+                0.9859501024515502,
+                -0.11222667947664411,
+                0.12372375636644917
+            ],
+            pixel_spacing=[1.0, 1.0],
+            spacing_between_slices=1.0,
+            coordinate_system='PATIENT'
+        ),
+        Volume.from_attributes(
+            array=np.random.randint(0, 100, (16, 16, 16)),
+            image_position=[0.0, 0.0, 0.0],
+            image_orientation=[
+                -0.4161787354389772,
+                0.440351232633788,
+                0.7955413578729377,
+                -0.2771687765466385,
+                -0.8947097563255662,
+                0.350245515665062
+            ],
+            pixel_spacing=[1.0, 1.0],
+            spacing_between_slices=1.0,
+            coordinate_system='PATIENT'
+        ),
+        # entirely random
+        Volume.from_attributes(
+            array=np.random.randint(113, 257, (192, 249, 84)),
+            image_position=[156.03935104, -57.61106994, -108.37601079],
+            image_orientation=[
+                -0.3056572521325831,
+                -0.9434667295206645,
+                -0.12823484123411946,
+                0.9440777770580763,
+                -0.3177984972478506,
+                0.08787073467366081
+            ],
+            pixel_spacing=[3.34201481, 2.35548103],
+            spacing_between_slices=2.82618053,
+            coordinate_system='PATIENT'
+        ),
+        Volume.from_attributes(
+            array=np.random.randint(180, 214, (96, 121, 50)),
+            image_position=[93.43769804, -184.44672839, -153.64700033],
+            image_orientation=[
+                -0.7034764816836532,
+                0.05445328347346446,
+                -0.7086294374614612,
+                -0.24703727011677554,
+                -0.9536256538386763,
+                0.1719613314498601
+            ],
+            pixel_spacing=[2.84370598, 0.69898499],
+            spacing_between_slices=0.57265037,
+            coordinate_system='PATIENT'
+        ),
+        Volume.from_attributes(
+            array=np.random.randint(81, 214, (34, 123, 59)),
+            image_position=[-252.23051789, 146.90528128, 84.40363858],
+            image_orientation=[
+                -0.7374895773326877,
+                0.5716640328411087,
+                0.35959610242811746,
+                -0.35709313761636013,
+                -0.7820074118766409,
+                0.510831575802926
+            ],
+            pixel_spacing=[3.52582689, 3.90364516],
+            spacing_between_slices=2.29457888,
+            coordinate_system='PATIENT'
+        ),
+
+    ]
+)
+def test_roundtrip(vol: Volume):
+    metatensor = vol.to_monai()
+
+    assert np.allclose(vol.get_affine('RAS'), metatensor.affine, atol=1e-4)
+    assert (vol.array == metatensor.numpy()).all()
+
+    monai_roundtrip = Volume.from_monai(metatensor)
+
+    assert np.allclose(vol.affine, monai_roundtrip.affine, atol=1e-4)
+    assert (vol.array == monai_roundtrip.array).all()
+
+
+@pytest.mark.parametrize(
+    'zip_url',
+    [
+        f'{DCM_QA_MPRAGE}/In/2_t1_mp2rage_sag_p3_32.zip',
+        f'{DCM_QA_MPRAGE}/In/5_HCP_T1.zip',
+    ]
+)
+def test_metatensor_equivalence_zip(zip_url: str):
+    vol, series, metatensor = read_github_zip_volume_and_metatensor(zip_url)
+
+    orientation = get_closest_patient_orientation(
+        convert_affine_to_convention(
+            metatensor.affine.numpy(),
+            from_reference_convention='RAS',
+            to_reference_convention='LPS'
+        )
+    )
+    oriented_vol = vol.to_patient_orientation(orientation)
+
+    assert np.allclose(
+        oriented_vol.get_affine('RAS'),
+        metatensor.affine.numpy(),
+        atol=1e-4
+    )
+    assert (oriented_vol.array == metatensor.numpy()).all()
+
+
+@pytest.mark.parametrize(
+    'dcm_urls',
+    [
+        [
+            f'{DCM_QA_ME}/In/2_me_FieldMap_GRE/{i:04d}.dcm'
+            for i in range(1, 37)
+        ],
+        [
+            f'{DCM_QA_PDT2}/In/Siemens/VE11/{i:04d}.dcm'
+            for i in range(1, 36)
+        ],
+    ]
+)
+def test_metatensor_equivalence_series(dcm_urls: Sequence[str]):
+    vol, series, metatensor = read_github_series_volume_and_metatensor(dcm_urls)
+
+    orientation = get_closest_patient_orientation(
+        convert_affine_to_convention(
+            metatensor.affine.numpy(),
+            from_reference_convention='RAS',
+            to_reference_convention='LPS'
+        )
+    )
+    oriented_vol = vol.to_patient_orientation(orientation)
+
+    assert np.allclose(
+        oriented_vol.get_affine('RAS'),
+        metatensor.affine.numpy(),
+        atol=1e-4
+    )
+    assert (oriented_vol.array == metatensor.numpy()).all()
+
+
+@pytest.mark.parametrize(
+    'segfile,channel_first,convert_to_ras,spacing,spatial_shape,affine,sum',
+    [
+        [
+            'seg_image_sm_control.dcm',
+            True,
+            True,
+            (1.0, 0.000499, 0.000499),
+            (1, 50, 50),
+            np.array([[0., 4.99e-4, 0., -23.449374],
+                      [0., 0., 4.99e-4, -25.691075],
+                      [1., 0., 0., 1.01],
+                      [0., 0., 0., 1.]]),
+            523
+        ],
+        [
+            'seg_image_sm_dots_tiled_full.dcm',
+            False,
+            True,
+            (1.0, 0.000499, 0.000499),
+            (1, 50, 50),
+            np.array([[0., 4.99e-4, 0., -23.449873],
+                      [0., 0., 4.99e-4, -25.691574],
+                      [1., 0., 0., 0.],
+                      [0., 0., 0., 1.]]),
+            200
+        ],
+        [
+            'seg_image_ct_true_fractional.dcm',
+            True,
+            False,
+            (1.25, 0.488281, 0.488281),
+            (3, 16, 16),
+            np.array([[0., 0., 0.488281, -125.],
+                      [0., 0.488281, 0., -128.100006],
+                      [-1.25, 0., 0., 105.519997],
+                      [0., 0., 0., 1.]]),
+            326.02353
+        ],
+        [
+            'seg_image_ct_binary_overlap.dcm',
+            False,
+            False,
+            (1.25, 0.488281, 0.488281),
+            (165, 16, 16),
+            np.array([[0., 0., 0.488281, -125.],
+                      [0., 0.488281, 0., -128.100006],
+                      [-1.25, 0., 0., 105.519997],
+                      [0., 0., 0., 1.]]),
+            80
+        ],
+        [
+            'seg_image_sm_numbers.dcm',
+            True,
+            True,
+            (1.0, 0.000499, 0.000499),
+            (1, 50, 50),
+            np.array([[0., 4.99e-4, 0., -23.449374],
+                      [0., 0., 4.99e-4, -25.691075],
+                      [1., 0., 0., 1.01],
+                      [0., 0., 0., 1.]]),
+            523
+        ],
+        [
+            'seg_image_ct_binary_fractional.dcm',
+            True,
+            True,
+            (1.25, 0.488281, 0.488281),
+            (3, 16, 16),
+            np.array([[0., 0., -0.488281, 125.],
+                      [0., -0.488281, 0., 128.100006],
+                      [-1.25, 0., 0., 105.519997],
+                      [0., 0., 0., 1.]]),
+            638.0
+        ],
+        [
+            'seg_image_ct_binary_single_frame.dcm',
+            True,
+            True,
+            (5.0, 0.661468, 0.661468),
+            (1, 128, 128),
+            np.array([[0., 0., -0.661468, 158.135803],
+                      [0., -0.661468, 0., 179.035797],
+                      [-5., 0., 0., -75.699997],
+                      [0., 0., 0., 1.]]),
+            1832
+        ],
+        [
+            'seg_image_sm_dots.dcm',
+            True,
+            True,
+            (1.0, 0.000499, 0.000499),
+            (1, 50, 50),
+            np.array([[0., 4.99e-4, 0., -23.449374],
+                      [0., 0., 4.99e-4, -25.691075],
+                      [1., 0., 0., 1.01],
+                      [0., 0., 0., 1.]]),
+            200
+        ],
+        [
+            'seg_image_sm_control_labelmap.dcm',
+            True,
+            True,
+            (1.0, 0.000499, 0.000499),
+            (1, 50, 50),
+            np.array([[0., 4.99e-4, 0., -23.449873],
+                      [0., 0., 4.99e-4, -25.691574],
+                      [1., 0., 0., 0.],
+                      [0., 0., 0., 1.]]),
+            523
+        ],
+        [
+            'seg_image_sm_control_labelmap_palette_color.dcm',
+            True,
+            True,
+            (1.0, 0.000499, 0.000499),
+            (1, 50, 50),
+            np.array([[0., 4.99e-4, 0., -23.449873],
+                      [0., 0., 4.99e-4, -25.691574],
+                      [1., 0., 0., 0.],
+                      [0., 0., 0., 1.]]),
+            523
+        ],
+        [
+            'seg_image_ct_binary.dcm',
+            True,
+            True,
+            (1.25, 0.488281, 0.488281),
+            (3, 16, 16),
+            np.array([[0., 0., -0.488281, 125.],
+                      [0., -0.488281, 0., 128.100006],
+                      [-1.25, 0., 0., 105.519997],
+                      [0., 0., 0., 1.]]),
+            638
+        ]
+    ]
+)
+def test_segmentation(
+    segfile,
+    channel_first,
+    convert_to_ras,
+    spacing,
+    spatial_shape,
+    affine,
+    sum
+):
+    seg = segread(TEST_DATA / segfile)
+    vol = seg.get_volume()
+    metatensor = vol.to_monai(
+        convert_to_ras=convert_to_ras,
+        ensure_channel_first=channel_first,
+    )
+    meta = metatensor.meta
+
+    if convert_to_ras:
+        space = monai.utils.enums.SpaceKeys.RAS
+
+    else:
+        space = monai.utils.enums.SpaceKeys.LPS
+
+    assert meta[MetaKeys.SPACE] == space
+    assert vol.spacing == meta[ImageStatsKeys.SPACING] == spacing
+    assert (np.array(vol.spatial_shape) == meta[MetaKeys.SPATIAL_SHAPE]).all()
+    assert (
+        metatensor.shape[1:] if channel_first else metatensor.shape[:3] ==
+        spatial_shape
+    )
+    assert np.allclose(vol.get_affine(space), affine, atol=1e-4)
+    assert np.allclose(metatensor.affine, affine, atol=1e-4)
+    assert np.allclose(meta[MetaKeys.ORIGINAL_AFFINE], affine, atol=1e-4)
+    assert np.allclose(meta[MetaKeys.AFFINE], affine, atol=1e-4)
+    assert meta[MetaKeys.ORIGINAL_CHANNEL_DIM] == -1
+    assert vol.array.sum() == metatensor.numpy().sum() == sum
+
+
+def test_channels():
+    array = np.zeros((10, 10, 10, 1))
+    volume = Volume.from_attributes(
+        array=array,
+        image_position=(0.0, 0.0, 0.0),
+        image_orientation=(1.0, 0.0, 0.0, 0.0, 1.0, 0.0),
+        pixel_spacing=(1.0, 1.0),
+        spacing_between_slices=2.0,
+        channels={'OpticalPathIdentifier': ['path1']},
+        coordinate_system="PATIENT",
+    )
+
+    metatensor = volume.to_monai(ensure_channel_first=False)
+    assert metatensor.shape == (10, 10, 10, 1)
+
+    array = np.zeros((10, 10, 10, 1))
+    volume = Volume.from_attributes(
+        array=array,
+        image_position=(0.0, 0.0, 0.0),
+        image_orientation=(1.0, 0.0, 0.0, 0.0, 1.0, 0.0),
+        pixel_spacing=(1.0, 1.0),
+        spacing_between_slices=2.0,
+        channels={'OpticalPathIdentifier': ['path1']},
+        coordinate_system="PATIENT",
+    )
+
+    metatensor = volume.to_monai(ensure_channel_first=True)
+    assert metatensor.shape == (1, 10, 10, 10)
+
+    array = np.zeros((10, 10, 10))
+    volume = Volume.from_attributes(
+        array=array,
+        image_position=(0.0, 0.0, 0.0),
+        image_orientation=(1.0, 0.0, 0.0, 0.0, 1.0, 0.0),
+        pixel_spacing=(1.0, 1.0),
+        spacing_between_slices=2.0,
+        coordinate_system="PATIENT",
+    )
+
+    metatensor = volume.to_monai(ensure_channel_first=False)
+    assert metatensor.shape == (10, 10, 10)
+
+    array = np.zeros((10, 10, 10))
+    volume = Volume.from_attributes(
+        array=array,
+        image_position=(0.0, 0.0, 0.0),
+        image_orientation=(1.0, 0.0, 0.0, 0.0, 1.0, 0.0),
+        pixel_spacing=(1.0, 1.0),
+        spacing_between_slices=2.0,
+        coordinate_system="PATIENT",
+    )
+
+    metatensor = volume.to_monai(ensure_channel_first=True)
+    assert metatensor.shape == (1, 10, 10, 10)
+
+    array = np.zeros((10, 10, 10, 2))
+    volume = Volume.from_attributes(
+        array=array,
+        image_position=(0.0, 0.0, 0.0),
+        image_orientation=(1.0, 0.0, 0.0, 0.0, 1.0, 0.0),
+        pixel_spacing=(1.0, 1.0),
+        spacing_between_slices=2.0,
+        channels={'OpticalPathIdentifier': ['path1', 'path2']},
+        coordinate_system="PATIENT",
+    )
+
+    volume.to_monai()
+
+    array = np.zeros((10, 10, 10, 1, 1))
+    volume = Volume.from_attributes(
+        array=array,
+        image_position=(0.0, 0.0, 0.0),
+        image_orientation=(1.0, 0.0, 0.0, 0.0, 1.0, 0.0),
+        pixel_spacing=(1.0, 1.0),
+        spacing_between_slices=2.0,
+        channels={
+            ChannelDescriptor('Channel0', True, str): ['class0'],
+            ChannelDescriptor('Channel1', True, str): ['class1']
+        },
+        coordinate_system="PATIENT",
+    )
+
+    with pytest.raises(
+        ValueError,
+        match=(
+            'Monai conversion does not currently support'
+            ' volumes with multiple channel dimensions.'
+        )
+    ):
+        volume.to_monai()
+
+    array = np.zeros((1, 10, 10, 10))
+    metatensor = monai.data.MetaTensor(array)
+    Volume.from_monai(metatensor)
+
+    array = np.zeros((10, 10, 10, 1))
+    metatensor = monai.data.MetaTensor(array)
+    Volume.from_monai(metatensor, channel_dim=-1)
+
+    array = np.zeros((1, 10, 10, 10))
+    metatensor = monai.data.MetaTensor(array)
+    Volume.from_monai(
+        metatensor,
+        channels={'OpticalPathIdentifier': ['path1']}
+    )
+
+    array = np.zeros((2, 10, 10, 10))
+    metatensor = monai.data.MetaTensor(array)
+    Volume.from_monai(
+        metatensor,
+        channels={'OpticalPathIdentifier': ['path1', 'path2']}
+    )
+
+    array = np.zeros((1, 1, 10, 10, 10))
+    metatensor = monai.data.MetaTensor(array)
+
+    with pytest.raises(
+        ValueError,
+        match=(
+            'Monai conversion does not currently support'
+            ' volumes with multiple channel dimensions.'
+        )
+    ):
+        Volume.from_monai(metatensor)
+
+    array = np.zeros((2, 10, 10, 10))
+    metatensor = monai.data.MetaTensor(array)
+
+    with pytest.raises(
+        ValueError,
+        match=(
+            'Monai conversion requires `channels` be specified'
+            ' for volumes with >=2 channels.'
+        )
+    ):
+        Volume.from_monai(metatensor)
+
+
+@pytest.mark.parametrize(
+    'volume,squeeze_dim',
+    [
+        [
+            Volume(
+                array=np.random.rand(1, 10, 10),
+                affine=np.array(
+                    [[0., 4.99e-4, 0., -23.449374],
+                     [0., 0., 4.99e-4, -25.691075],
+                     [1., 0., 0., 1.01],
+                     [0., 0., 0., 1.]]
+                ),
+                coordinate_system="PATIENT",
+            ),
+            0
+        ],
+        [
+            Volume(
+                array=np.random.rand(10, 1, 10),
+                affine=np.array(
+                    [[0., 4.99e-04, 0., -23.449873],
+                     [0., 0., 4.99e-04, -25.691574],
+                     [1., 0., 0., 0.],
+                     [0., 0., 0., 1.]]
+                ),
+                coordinate_system="PATIENT",
+            ),
+            1
+        ],
+        [
+            Volume(
+                array=np.random.rand(10, 10, 1),
+                affine=np.array(
+                    [[0., 4.99e-04, 0., -23.449873],
+                     [0., 0., 4.99e-04, -25.691574],
+                     [1., 0., 0., 0.],
+                     [0., 0., 0., 1.]]
+                ),
+                coordinate_system="PATIENT",
+            ),
+            2
+        ],
+        [
+            Volume(
+                array=np.random.rand(1, 10, 10),
+                affine=np.array(
+                    [[0.51484336, -2.45034467, -0.29838402, 31.96125748],
+                     [0.84430051, 1.01480423, -1.0149242, 19.85095778],
+                     [1.16995367, 0.34594871, 0.86372826, -39.22734667],
+                     [0., 0., 0., 1.]]
+                ),
+                coordinate_system="PATIENT",
+            ),
+            0
+        ],
+        [
+            Volume(
+                array=np.random.rand(10, 1, 10),
+                affine=np.array(
+                    [[-6.62073016e-03, -0.648063598, 6.39172474e-02, -25.55167],
+                     [-0.786401952, -1.4152762e-02, -1.50316788, -12.367773],
+                     [0.950988641, -1.621514e-02, -1.24257108, -11.5527702],
+                     [0., 0., 0., 1.]]
+                ),
+                coordinate_system="PATIENT",
+            ),
+            1
+        ],
+        [
+            Volume(
+                array=np.random.rand(10, 10, 1),
+                affine=np.array(
+                    [[0.783406238, 0.312050746, 2.29152341e-02, 45.850208],
+                     [7.69377955e-02, -0.998526806, -0.712275579, -49.9411473],
+                     [-9.27305439e-02, 1.80779755, -0.397376895, 18.0307572],
+                     [0., 0., 0., 1.]]
+                ),
+                coordinate_system="PATIENT",
+            ),
+            2
+        ],
+        [
+            Volume(
+                array=np.random.rand(1, 10, 10),
+                affine=np.array(
+                    [[-1.8249537, 0.357257741, -0.365592144, 35.5553115],
+                     [1.06765446, 0.586750214, -0.637285403, 18.5404252],
+                     [-2.58571459e-02, -0.987438965, -0.510956107, 27.5584146],
+                     [0., 0., 0., 1.]]
+                ),
+                coordinate_system="PATIENT",
+            ),
+            0
+        ],
+        [
+            Volume(
+                array=np.random.rand(10, 1, 10),
+                affine=np.array(
+                    [[0.40659061, -0.14035702, 2.52546278, -22.24126432],
+                     [-2.01663853, -0.32789655, 0.12001713, -27.24361216],
+                     [1.09727089, -0.55062153, -0.71522747, -25.00163784],
+                     [0., 0., 0., 1.]]
+                ),
+                coordinate_system="PATIENT",
+            ),
+            1
+        ],
+        [
+            Volume(
+                array=np.random.rand(10, 10, 1),
+                affine=np.array(
+                    [[-0.692324440, 1.81254936, 0.770028884, 19.9835788],
+                     [-1.01187701, -1.29172094, 1.06442359, 4.6090612],
+                     [2.01112784, -2.59508652e-02, 0.800633129, -38.2126518],
+                     [0., 0., 0., 1.]]
+                ),
+                coordinate_system="PATIENT",
+            ),
+            2
+        ],
+    ]
+)
+def test_squeeze_affine(volume, squeeze_dim):
+    metatensor = volume.to_monai(squeeze_dim=squeeze_dim)
+
+    assert (
+        volume.array.squeeze(squeeze_dim) == metatensor.numpy()
+    ).all()
+    validate_affine(
+        affine_3d=volume.affine,
+        affine_2d=metatensor.affine.numpy(),
+        squeeze_dim=squeeze_dim,
+        shape=volume.shape
+    )
+
+
+def test_squeeze_errors():
+    volume = Volume(
+        array=np.random.rand(1, 10, 10),
+        affine=np.array([[0., 4.99e-4, 0., -23.449374],
+                         [0., 0., 4.99e-4, -25.691075],
+                         [1., 0., 0., 1.01],
+                         [0., 0., 0., 1.]]),
+        coordinate_system="PATIENT",
+    )
+    squeeze_dim = 3
+
+    with pytest.raises(
+        ValueError,
+        match=re.escape(
+            'If provided, `squeeze_dim` must be a spatial'
+            ' dimension (0, 1, 2).'
+        )
+    ):
+        volume.to_monai(squeeze_dim=squeeze_dim)
+
+    volume = Volume(
+        array=np.random.rand(10, 10, 10),
+        affine=np.array([[0., 4.99e-4, 0., -23.449374],
+                         [0., 0., 4.99e-4, -25.691075],
+                         [1., 0., 0., 1.01],
+                         [0., 0., 0., 1.]]),
+        coordinate_system="PATIENT",
+    )
+    squeeze_dim = 0
+
+    with pytest.raises(
+        ValueError,
+        match=re.escape(
+            f'`squeeze_dim={squeeze_dim}` does not correspond'
+            ' to a singleton dimension. Array has shape:'
+            f' {volume.array.shape}.'
+        )
+    ):
+        volume.to_monai(squeeze_dim=squeeze_dim)
diff --git a/tests/utils.py b/tests/utils.py
index 31a4f2c0..dfdb55c5 100644
--- a/tests/utils.py
+++ b/tests/utils.py
@@ -19,6 +19,7 @@
 DCM_QA_MPRAGE = 'https://github.com/neurolabusc/dcm_qa_mprage/raw/refs/heads/main'  # noqa: E501
 DCM_QA_ME = 'https://github.com/neurolabusc/dcm_qa_me/raw/refs/heads/master'
 DCM_QA_PDT2 = 'https://github.com/neurolabusc/dcm_qa_pdt2/raw/refs/heads/main'
+TEST_DATA = Path(__file__).parent.parent / "data" / "test_files"
 
 
 def write_and_read_dataset(dataset: Dataset):