Functions
Tools for aligning spatial transcriptomics data
Functions:
|
Run LDDMM between a pair of images. |
|
LDDMM for 3D to 2D slice mapping. |
|
Compute an affine transformation from points. |
|
Create dataframe with region annotations and 3D coordinates of target brain slice |
|
Create sample points to transform source to target from affine and velocity. |
|
Create sample points to transform source to target from affine and velocity. |
|
clip an arrays values between 0 and 1. |
|
Create 3D altas image and region annotations |
|
Create 3D altas ontology |
|
Given a set of pixel locations, returns an extent 4-tuple for use with np.imshow. |
|
Interpolate the 2D image I, with regular grid positions stored in x (1d arrays), at the positions stored in phii (3D arrays with first channel storing component) |
|
Interpolate the 3D image I, with regular grid positions stored in x (1d arrays), at the positions stored in phii (4D arrays with first channel storing component) |
|
Create a scree plot for a given set of eigenvalues. |
|
Linearly normalizes an array between two specifed values. |
Plot all brain regions in target brain slice with different color. |
|
|
Plot subset of brain regions in target brain slice with different color. |
|
Rasterize a spatial transcriptomics dataset into a density image |
|
Rasterize a spatial transcriptomics dataset into a density image and perform PCA on it. |
|
Rasterize a spatial transcriptomics dataset into a density image |
|
Create an RGB image principal components 1-3 and 4-6. |
|
Convert a linear transform matrix and a translation vector into an affine matrix. |
|
Convert a linear transform matrix and a translation vector into an affine matrix. |
|
Transform an image |
|
Transform an image with an affine matrix |
|
Transform an image |
|
Transform points. |
|
Transform points. |
|
Integrate a velocity field over time to return a position field (diffeomorphism). |
|
Integrate a 3D velocity field over time to return a 3D position field (diffeomorphism). |
- STalign.STalign.LDDMM(xI, I, xJ, J, pointsI=None, pointsJ=None, L=None, T=None, A=None, v=None, xv=None, a=500.0, p=2.0, expand=2.0, nt=3, niter=5000, diffeo_start=0, epL=2e-08, epT=0.2, epV=2000.0, sigmaM=1.0, sigmaB=2.0, sigmaA=5.0, sigmaR=500000.0, sigmaP=20.0, device='cpu', dtype=torch.float64, muB=None, muA=None)
Run LDDMM between a pair of images.
This jointly estimates an affine transform A, and a diffeomorphism phi. The map is off the form x -> A phi x
Parameters
- xIlist of torch tensor
Location of voxels in source image I
- Itorch tensor
source image I, with channels along first axis
- xJlist of torch tensor
Location of voxels in target image J
- Jtorch tensor
Target image J, with channels along first axis
- Ltorch tensor
Initial guess for linear transform (2x2 torch tensor). Defaults to None (identity).
- Ttorch tensor
Initial guess for translation (2 element torch tensor). Defaults to None (identity)
- Atorch tensor
Initial guess for affine matrix. Either L and T can be specified, or A, but not both. Defaults to None (identity).
- vtorch tensor
Initial guess for velocity field
- xvtorch tensor
pixel locations for velocity field
- afloat
Smoothness scale of velocity field (default 500.0)
- pfloat
Power of Laplacian in velocity regularization (default 2.0)
- expandfloat
Factor to expand size of velocity field around image boundaries (default 2.0)
- ntint
Number of timesteps for integrating velocity field (default 3). Ignored if you input v.
- pointsItorch tensor
N x 2 set of corresponding points for matching in source image. Default None (no points).
- pointsJtorch tensor
N x 2 set of corresponding points for matching in target image. Default None (no points).
- niterint
Number of iterations of gradient descent optimization
- diffeo_startint
Number of iterations of gradient descent optimization for affine only, before nonlinear deformation.
- epLfloat
Gradient descent step size for linear part of affine.
- epTfloat
Gradient descent step size of translation part of affine.
- epVfloat
Gradient descent step size for velocity field.
- sigmaMfloat
Standard deviation of image matching term for Gaussian mixture modeling in cost function. This term generally controls matching accuracy with smaller corresponding to more accurate. As an common example (rule of thumb), you could chose this parameter to be the variance of the pixels in your target image.
- sigmaBfloat
Standard deviation of backtround term for Gaussian mixture modeling in cost function. If there is missing tissue in target, we may label some pixels in target as background, and not enforce matching here.
- sigmaAfloat
Standard deviation of artifact term for Gaussian mixture modeling in cost function. If there are artifacts in target or other lack of corresponding between template and target, we may label some pixels in target as artifact, and not enforce matching here.
- sigmaR: float
Standard deviation for regularization. Smaller sigmaR means a smoother resulting transformation. Regularization is of the form: 0.5/sigmaR^2 int_0^1 int_X |Lv|^2 dx dt.
- sigmaP: float
Standard deviation for matching of points. Cost is of the form 0.5/sigmaP^2 sum_i (source_point_i - target_point_i)^2
- device: str
torch device. defaults to ‘cpu’. Can also be ‘cuda:0’ for example.
- dtype: torch dtype
torch data type. defaults to torch.float64
- muA: torch tensor whose dimension is the same as the target image
Defaults to None, which means we estimate this. If you provide a value, we will not estimate it. If the target is a RGB image, this should be a tensor of size 3. If the target is a grayscale image, this should be a tensor of size 1.
- muB: torch tensor whose dimension is the same as the target image
Defaults to None, which means we estimate this. If you provide a value, we will not estimate it.
Returns a dictionary
{ ‘A’: torch tensor
Affine transform
- ‘v’: torch tensor
Velocity field
- ‘xv’: list of torch tensor
Pixel locations in v
- ‘WM’: torch tensor
Resulting weight 2D array (matching)
- ‘WB’: torch tensor
Resulting weight 2D array (background)
- ‘WA’: torch tensor
Resulting weight 2D array (artifact)
}
- STalign.STalign.LDDMM_3D_to_slice(xI, I, xJ, J, pointsI=None, pointsJ=None, L=None, T=None, A=None, v=None, xv=None, a=500.0, p=2.0, expand=1.25, nt=3, niter=5000, diffeo_start=0, epL=1e-06, epT=10.0, epV=1000.0, sigmaM=1.0, sigmaB=2.0, sigmaA=5.0, sigmaR=100000000.0, sigmaP=20.0, device='cpu', dtype=torch.float64, muA=None, muB=None)
LDDMM for 3D to 2D slice mapping.
- muA: torch tensor whose dimension is the same as the target image
Defaults to None, which means we estimate this. If you provide a value, we will not estimate it. If the target is a RGB image, this should be a tensor of size 3. If the target is a grayscale image, this should be a tensor of size 1.
- muB: torch tensor whose dimension is the same as the target image
Defaults to None, which means we estimate this. If you provide a value, we will not estimate it.
- STalign.STalign.L_T_from_points(pointsI, pointsJ)
Compute an affine transformation from points.
Note for an affine transformation (6dof) we need 3 points.
Outputs, L,T should be rconstructed blockwize like [L,T;0,0,1]
Parameters
- pointsIarray
An Nx2 array of floating point numbers describing source points in ROW COL order (not xy)
- pointsJarray
An Nx2 array of floating point numbers describing target points in ROW COL order (not xy)
Returns
- Larray
A 2x2 linear transform array.
- Tarray
A 2 element translation vector
- STalign.STalign.analyze3Dalign(labelfile, xv, v, mat, xJ, dx, scale_x, scale_y, x, y, X_, Y_, namesdict, device='cpu')
Create dataframe with region annotations and 3D coordinates of target brain slice
- labelfilearray
File name of 3D atlas with region annotations.
- xvlist of array
Sample points for velocity
- varray
time dependent velocity field
- Aarray
Affine transformation matrix
- directionchar
‘f’ for forward and ‘b’ for backward. ‘b’ is default and is used for transforming images. ‘f’ is used for transforming points.
- XJarray
Sample points for target (meshgrid with ij index style). Defaults to None
to keep sampling on the xv.
- dxint
Step of rasterized image.
- scale_xdouble
Value that x positions of atlas were scaled by
- scale_ydouble
Value that y positions of atlas were scaled by
- xarray
X positions of target brain slice
- yarray
Y positions of target brain slice
- X_array
- Y_array
- namesdictDictionary
Dictionary of brain regions corresponding to region numbers
df : Dataframe containing each cell in original target brain slice, region annotations and 3D coordinates in terms of the altas.
- STalign.STalign.build_transform(xv, v, A, direction='b', XJ=None)
Create sample points to transform source to target from affine and velocity.
Parameters
- xvlist of array
Sample points for velocity
- varray
time dependent velocity field
- Aarray
Affine transformation matrix
- directionchar
‘f’ for forward and ‘b’ for backward. ‘b’ is default and is used for transforming images. ‘f’ is used for transforming points.
- XJarray
Sample points for target (meshgrid with ij index style). Defaults to None to keep sampling on the xv.
Returns
- Xsarray
Sample points in mehsgrid format.
- STalign.STalign.build_transform3D(xv, v, A, direction='b', XJ=None)
Create sample points to transform source to target from affine and velocity.
Parameters
- xvlist of array
Sample points for velocity
- varray
time dependent velocity field
- Aarray
Affine transformation matrix
- directionchar
‘f’ for forward and ‘b’ for backward. ‘b’ is default and is used for transforming images. ‘f’ is used for transforming points.
- XJarray
Sample points for target (meshgrid with ij index style). Defaults to None to keep sampling on the xv.
Returns
- Xsarray
Sample points in mehsgrid format.
- STalign.STalign.clip(I)
clip an arrays values between 0 and 1. Useful for visualizing
Parameters
- Itorch tensor
A torch tensor to clip.
Returns
- Ictorch tensor
Clipped torch tensor.
- STalign.STalign.download_aba_image_labels(imageurl, labelurl, imagefile, labelfile)
Create 3D altas image and region annotations
Parameters
- imageurlarray
Link containing atlas cell stained image.
- labelurlarray
Link containing atlas altas images for each voxel in imageurl.
- imagefilearray
File location to save imageurl information.
- labelfilearray
File location to save labelurl information.
Returns
- imagefilearray
File location to save imageurl information.
- labelfilearray
File location to save labelurl information.
- STalign.STalign.download_aba_ontology(url, file_name)
Create 3D altas ontology
Parameters
url : Link to url contain atlas ontology file_name : File name to save atlas ontology
Returns
ontology_name : File name to store ontologies namesdict : Dictionary of all brain rgeion names
- STalign.STalign.extent_from_x(xJ)
Given a set of pixel locations, returns an extent 4-tuple for use with np.imshow.
Note inputs are locations of pixels along each axis, i.e. row column not xy.
Parameters
- xJlist of torch tensors
Location of pixels along each axis
Returns
- extenttuple
(xmin, xmax, ymin, ymax) tuple
Examples
>>> extent_from_x(xJ) >>> fig,ax = plt.subplots() >>> ax.imshow(J,extent=extentJ)
- STalign.STalign.interp(x, I, phii, **kwargs)
Interpolate the 2D image I, with regular grid positions stored in x (1d arrays), at the positions stored in phii (3D arrays with first channel storing component)
Parameters
- xlist of arrays
List of arrays storing the pixel locations along each image axis. convention is row column order not xy.
- Iarray
Image array. First axis should contain different channels.
- phiiarray
Sampling array. First axis should contain sample locations corresponding to each axis.
- **kwargsdict
Other arguments fed into the torch interpolation function torch.nn.grid_sample
Returns
- outtorch tensor
The image I resampled on the points defined in phii.
Notes
Convention is to use align_corners=True.
This uses the torch library.
- STalign.STalign.interp3D(x, I, phii, **kwargs)
Interpolate the 3D image I, with regular grid positions stored in x (1d arrays), at the positions stored in phii (4D arrays with first channel storing component)
Parameters
- xlist of arrays
List of arrays storing the pixel locations along each image axis. convention is row column order not xy.
- Iarray
Image array. First axis should contain different channels.
- phiiarray
Sampling array. First axis should contain sample locations corresponding to each axis.
- **kwargsdict
Other arguments fed into the torch interpolation function torch.nn.grid_sample
Returns
- outtorch tensor
The image I resampled on the points defined in phii.
Notes
Convention is to use align_corners=True.
This uses the torch library.
- STalign.STalign.make_scree(W, name, p=6)
Create a scree plot for a given set of eigenvalues.
Parameters
- Wnumpy array
Eigenvalues in descending order
- namestr
Name of the dataset the eigenvalue originate from Will be included in plot title
- pint
Number of eigenvalue from W to plot. Defaults to 6
Returns
- figmatplotlib Figure
Figure handle for scree plot
Raises
- ValueError
If p is larger than the number to eigenvalues in W.
- STalign.STalign.normalize(arr, t_min=0, t_max=1)
Linearly normalizes an array between two specifed values.
Parameters
- arrnumpy array
array to be normalized
- t_minint or float
Lower bound of normalization range
- t_maxint or float
Upper bound of normalization range
Returns
- norm_arrnumpy array
1D array with normalized arr values
- STalign.STalign.plot_brain_regions(df)
Plot all brain regions in target brain slice with different color.
Parameters
- dfDataFrame
Dataframe containing each cell in original target brain slice, region annotations and 3D coordinates in terms of the altas.
Returns
None
- STalign.STalign.plot_subset_brain_regions(df, brain_regions)
Plot subset of brain regions in target brain slice with different color.
Parameters
- dfDataFrame
Dataframe containing each cell in original target brain slice, region annotations and 3D coordinates in terms of the altas.
- brain_regionsarray
Subset of brain regions of interest (i.e. [‘CA1’, ‘CP’, ‘DG-sg’])
Returns
None
- STalign.STalign.rasterize(x, y, g=array([1.]), dx=30.0, blur=1.0, expand=1.1, draw=10000, wavelet_magnitude=False, use_windowing=True)
Rasterize a spatial transcriptomics dataset into a density image
Paramters
- xnumpy array of length N
x location of cells
- ynumpy array of length N
y location of cells
- gnumpy array of length N
RNA count of cells If not given, density image is created
- dxfloat
Pixel size to rasterize data (default 30.0, in same units as x and y)
- blurfloat or list of floats
Standard deviation of Gaussian interpolation kernel. Units are in number of pixels. Can be aUse a list to do multi scale.
- expandfloat
Factor to expand sampled area beyond cells. Defaults to 1.1.
- drawint
If True, draw a figure every draw points return its handle. Defaults to False (0).
- wavelet_magnitudebool
If True, take the absolute value of difference between scales for raster images. When using this option blur should be sorted from greatest to least.
Returns
- Xnumpy array
Locations of pixels along the x axis
- Ynumpy array
Locations of pixels along the y axis
- Mnumpy array
A rasterized image with len(blur) channels along the first axis
- figmatplotlib figure handle
If draw=True, returns a figure handle to the drawn figure.
Raises
- Exception
If wavelet_magnitude is set to true but blur is not sorted from greatest to least.
Examples
Rasterize a dataset at 30 micron pixel size, with three kernels.
>>> X,Y,M,fig = tools.rasterize(x,y,dx=30.0,blur=[2.0,1.0,0.5],draw=10000)
Rasterize a dataset at 30 micron pixel size, with three kernels, using difference between scales.
>>> X,Y,M,fig = tools.rasterize(x,y,dx=30.0,blur=[2.0,1.0,0.5],draw=10000, wavelet_magnitude=True)
- STalign.STalign.rasterizePCA(x, y, G)
Rasterize a spatial transcriptomics dataset into a density image and perform PCA on it.
Parameters
- xnumpy array of length N
x location of cells
- ynumpy array of length N
y location of cells
- Gpandas Dataframe of shape (N,M)
gene expression level of cells for M genes
Returns
- Xnumpy array
A rasterized image for each gene with the channel along the first axis
- Ynumpy array
A rasterized image unraveled to 1-D for each gene; data is centered
- Wnumpy array
The eigenvalues for the principal components in descending order
- Vnumpy array
The normalized eigenvectors for the principal components V[:, i] corresponds to eigenvalue at W[i]
- Znumpy array
X rotated to align with the principal component axes
- nrowsint
Row dimension of rasterized image
- ncolsint
Column dimension of rasterized image
Notes
Each value/row at the same index in x, y, and G should all correspond to the same cell. x[i] <-> y[i] <-> G[i,:]
- STalign.STalign.rasterize_with_signal(x, y, s=None, dx=30.0, blur=1.0, expand=1.1, draw=0, wavelet_magnitude=False, use_windowing=True)
Rasterize a spatial transcriptomics dataset into a density image
Paramters
- xnumpy array of length N
x location of cells
- ynumpy array of length N
y location of cells
- snumpy array of length N by M for M signals
signal value should be length NxM
- dxfloat
Pixel size to rasterize data (default 30.0, in same units as x and y)
- blurfloat or list of floats
Standard deviation of Gaussian interpolation kernel. Units are in number of pixels. Can be aUse a list to do multi scale.
- expandfloat
Factor to expand sampled area beyond cells. Defaults to 1.1.
- drawint
If True, draw a figure every draw points return its handle. Defaults to False (0).
- wavelet_magnitudebool
If True, take the absolute value of difference between scales for raster images. When using this option blur should be sorted from greatest to least.
Returns
- Xnumpy array
Locations of pixels along the x axis
- Ynumpy array
Locations of pixels along the y axis
- Mnumpy array
A rasterized image with len(blur) channels along the last axis
- figmatplotlib figure handle
If draw=True, returns a figure handle to the drawn figure.
Raises
- Exception
If wavelet_magnitude is set to true but blur is not sorted from greatest to least.
Examples
Rasterize a dataset at 30 micron pixel size, with three kernels.
>>> X,Y,M,fig = tools.rasterize(x,y,dx=30.0,blur=[2.0,1.0,0.5],draw=10000)
Rasterize a dataset at 30 micron pixel size, with three kernels, using difference between scales.
>>> X,Y,M,fig = tools.rasterize(x,y,dx=30.0,blur=[2.0,1.0,0.5],draw=10000, wavelet_magnitude=True)
- STalign.STalign.saveRasters(X, Y, W, V, Z, scree_fig, nrows, ncols, name, path)
Create an RGB image principal components 1-3 and 4-6. Save numpy arrays necessary for rasterizing and performing PCA. Save scree plot and RGB image as jpeg, and their figure handles.
Parameters
- Xnumpy array
A rasterized image for each gene with the channel along the first axis
- Ynumpy array
A rasterized image unraveled to 1-D for each gene; data is centered
- Wnumpy array
The eigenvalues for the principal components in descending order
- Vnumpy array
The normalized eigenvectors for the principal components V[:, i] corresponds to eigenvalue at W[i]
- Znumpy array
X rotated to align with the principal component axes
- scree_figmatplotlib Figure
Figure handle for scree plot
- nrowsint
Row dimension of rasterized image
- ncolsint
Column dimension of rasterized image
- namestr
Name of the dataset the arrays originate from
- pathstr
Absolute path to folder where files should be saved
Returns
None
- STalign.STalign.to_A(L, T)
Convert a linear transform matrix and a translation vector into an affine matrix.
Parameters
- Ltorch tensor
2x2 linear transform matrix
- Ttorch tensor
2 element translation vector (note NOT 2x1)
Returns
- Atorch tensor
Affine transform matrix
- STalign.STalign.to_A_3D(L, T)
Convert a linear transform matrix and a translation vector into an affine matrix.
Parameters
- Ltorch tensor
3x3 linear transform matrix
- Ttorch tensor
3 element translation vector (note NOT 2x1)
Returns
- Atorch tensor
Affine transform matrix
- STalign.STalign.transform_image_source_to_target(xv, v, A, xI, I, XJ=None)
Transform an image
- STalign.STalign.transform_image_source_with_A(A, XI, I, XJ)
Transform an image with an affine matrix
Parameters
- Atorch tensor
Affine transform matrix
- XIlist of numpy arrays
List of arrays storing the pixel location in image I along each image axis. convention is row column order not xy. i.e, locations of pixels along the y axis (rows) followed by locations of pixels along the x axis (columns)
- Inumpy array
A rasterized image with len(blur) channels along the first axis
- XJlist of numpy arrays
List of arrays storing the pixel location in image I along each image axis. convention is row column order not xy. i.e, locations of pixels along the y axis (rows) followed by locations of pixels along the x axis (columns)
Returns
- AItorch tensor
image I after affine transformation A, with channels along first axis
- STalign.STalign.transform_image_target_to_source(xv, v, A, xJ, J, XI=None)
Transform an image
- STalign.STalign.transform_points_source_to_target(xv, v, A, pointsI)
Transform points. Note points are in row column order, not xy.
- STalign.STalign.transform_points_target_to_source(xv, v, A, pointsI)
Transform points. Note points are in row column order, not xy.
- STalign.STalign.v_to_phii(xv, v)
Integrate a velocity field over time to return a position field (diffeomorphism).
Parameters
- xvlist of torch tensor
List of 1D tensors describing locations of sample points in v
- vtorch tensor
5D (nt,2,v0,v1) velocity field
Returns
- phii: torch tensor
Inverse map (position field) computed by method of characteristics
- STalign.STalign.v_to_phii_3D(xv, v)
Integrate a 3D velocity field over time to return a 3D position field (diffeomorphism).
Parameters
- xvlist of torch tensor
List of 1D tensors describing locations of sample points in v
- vtorch tensor
5D (nt,3,v0,v1,v2) velocity field
Returns
- phii: torch tensor
Inverse map (position field) computed by method of characteristics