Note

Go to the end to download the full example code

# Using an automated approach to coregistration#

This example shows how to use the coregistration functions to perform an automated MEG-MRI coregistration via scripting. Generally the results of this approach are consistent with those obtained from manual coregistration [1].

Warning

The quality of the coregistration depends heavily upon the quality of the head shape points (HSP) collected during subject prepration and the quality of your T1-weighted MRI. Use with caution and check the coregistration error.

```
# Author: Jon Houck <jon.houck@gmail.com>
# Guillaume Favelier <guillaume.favelier@gmail.com>
#
# License: BSD-3-Clause
import numpy as np
import mne
from mne.coreg import Coregistration
from mne.io import read_info
data_path = mne.datasets.sample.data_path()
# data_path and all paths built from it are pathlib.Path objects
subjects_dir = data_path / 'subjects'
subject = 'sample'
fname_raw = data_path / 'MEG' / subject / f'{subject}_audvis_raw.fif'
info = read_info(fname_raw)
plot_kwargs = dict(subject=subject, subjects_dir=subjects_dir,
surfaces="head-dense", dig=True, eeg=[],
meg='sensors', show_axes=True,
coord_frame='meg')
view_kwargs = dict(azimuth=45, elevation=90, distance=0.6,
focalpoint=(0., 0., 0.))
```

```
Read a total of 3 projection items:
PCA-v1 (1 x 102) idle
PCA-v2 (1 x 102) idle
PCA-v3 (1 x 102) idle
```

## Set up the coregistration model#

```
fiducials = "estimated" # get fiducials from fsaverage
coreg = Coregistration(info, subject, subjects_dir, fiducials=fiducials)
fig = mne.viz.plot_alignment(info, trans=coreg.trans, **plot_kwargs)
```

```
Triangle neighbors and vertex normals...
Using high resolution head model in /home/circleci/mne_data/MNE-sample-data/subjects/sample/surf/lh.seghead
Triangle neighbors and vertex normals...
Estimating fiducials from fsaverage.
Using lh.seghead for head surface.
Channel types:: grad: 203, mag: 102
```

## Initial fit with fiducials#

Do first a coregistration fit using only 3 fiducial points. This allows to find a good initial solution before further optimization using head shape points. This can also be useful to detect outlier head shape points which are too far from the skin surface. One can see for example that on this dataset there is one such point and we will omit it from the subsequent fit.

```
coreg.fit_fiducials(verbose=True)
fig = mne.viz.plot_alignment(info, trans=coreg.trans, **plot_kwargs)
```

```
Aligning using fiducials
Start median distance: 17.29 mm
End median distance: 11.11 mm
Using lh.seghead for head surface.
Channel types:: grad: 203, mag: 102
```

## Refining with ICP#

Next we refine the transformation using a few iteration of the Iterative Closest Point (ICP) algorithm. As the initial fiducials are obtained from fsaverage and not from precise manual picking in the GUI we do a fit with reduced weight for the nasion.

```
coreg.fit_icp(n_iterations=6, nasion_weight=2., verbose=True)
fig = mne.viz.plot_alignment(info, trans=coreg.trans, **plot_kwargs)
```

```
Aligning using ICP
Start median distance: 11.11 mm
ICP 1 median distance: 6.57 mm
ICP 2 median distance: 4.39 mm
ICP 3 median distance: 3.48 mm
ICP 4 median distance: 2.86 mm
ICP 5 median distance: 2.43 mm
ICP 6 median distance: 2.14 mm
End median distance: 2.14 mm
Using lh.seghead for head surface.
Channel types:: grad: 203, mag: 102
```

## Omitting bad points#

It is now very clear that we have one point that is an outlier and that should be removed.

```
coreg.omit_head_shape_points(distance=5. / 1000) # distance is in meters
```

```
Coregistration: Excluding 3 head shape points with distance >= 0.005 m.
```

## Final coregistration fit#

```
coreg.fit_icp(n_iterations=20, nasion_weight=10., verbose=True)
fig = mne.viz.plot_alignment(info, trans=coreg.trans, **plot_kwargs)
mne.viz.set_3d_view(fig, **view_kwargs)
dists = coreg.compute_dig_mri_distances() * 1e3 # in mm
print(
f"Distance between HSP and MRI (mean/min/max):\n{np.mean(dists):.2f} mm "
f"/ {np.min(dists):.2f} mm / {np.max(dists):.2f} mm"
)
```

```
Aligning using ICP
Start median distance: 2.00 mm
ICP 1 median distance: 1.99 mm
ICP 2 median distance: 2.02 mm
ICP 3 median distance: 1.79 mm
ICP 4 median distance: 1.78 mm
ICP 5 median distance: 1.77 mm
End median distance: 1.77 mm
Using lh.seghead for head surface.
Channel types:: grad: 203, mag: 102
Distance between HSP and MRI (mean/min/max):
2.00 mm / 0.19 mm / 6.37 mm
```

Warning

Don’t forget to save the resulting `trans`

matrix!

```
mne.write_trans('/path/to/filename-trans.fif', coreg.trans)
```

Note

The `mne.coreg.Coregistration`

class has the ability to
compute MRI scale factors using
`set_scale_mode()`

that is useful
for creating surrogate MRI subjects, i.e., using a template MRI
(such as one from `mne.datasets.fetch_infant_template()`

)
matched to a subject’s head digitization. When scaling is desired,
a scaled surrogate MRI should be created using
`mne.scale_mri()`

.

## References#

**Total running time of the script:** ( 0 minutes 29.171 seconds)

**Estimated memory usage:** 51 MB