-
Notifications
You must be signed in to change notification settings - Fork 0
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
Colocalization module #17
Merged
Changes from 8 commits
Commits
Show all changes
9 commits
Select commit
Hold shift + click to select a range
ac6bcc2
Colocalization module
tdrose 3968687
Documentation for coloc module
tdrose a12792f
Docs examples fix.
tdrose 1de5c51
renaming of coloc_ml_processing, testing, and moving colocs to varp
tdrose 7bc4dca
Changed test scenarios.
tdrose 94a2203
Testing median filter for all ions of one pixel.
tdrose 7561423
Separated quantile thresholding and median filtering in test.
tdrose 178e56c
Added layer argument to anndata_to_image_array function.
tdrose 9d09565
Fixed typo
tdrose File filter
Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
There are no files selected for viewing
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -1,5 +1,11 @@ | ||
from metaspace_converter import colocalization | ||
from metaspace_converter.anndata_to_array import anndata_to_image_array | ||
from metaspace_converter.to_anndata import metaspace_to_anndata | ||
from metaspace_converter.to_spatialdata import metaspace_to_spatialdata | ||
|
||
__all__ = ["metaspace_to_anndata", "metaspace_to_spatialdata", "anndata_to_image_array"] | ||
__all__ = [ | ||
"metaspace_to_anndata", | ||
"metaspace_to_spatialdata", | ||
"anndata_to_image_array", | ||
"colocalization", | ||
] |
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -0,0 +1,117 @@ | ||
from typing import Optional, Tuple | ||
|
||
import numpy as np | ||
from anndata import AnnData | ||
from scipy import ndimage | ||
|
||
from metaspace_converter.anndata_to_array import anndata_to_image_array | ||
from metaspace_converter.constants import COLOCALIZATION | ||
|
||
|
||
def coloc_ml_preprocessing( | ||
adata: AnnData, | ||
layer: Optional[str] = "coloc_ml_preprocessing", | ||
median_filter_size: Tuple[int, int] = (3, 3), | ||
quantile_threshold: float = 0.5, | ||
): | ||
""" | ||
Preprocessing for colocalization analysis according to the colocML publication | ||
(https://doi.org/10.1093/bioinformatics/btaa085). | ||
|
||
In the publication the authors evaluated colocalization metrics and preprocessing approaches. | ||
They found the best performance for | ||
1) median filtering of ion images with a (3, 3) kernel size and | ||
2) quantile thresholding ad 50%, meaning all pixels with intensities below the 50% | ||
quantile set to 0. | ||
|
||
This function performs the same preprocessing steps. | ||
Recommended for call before running the ``colocalization`` function. | ||
|
||
Args: | ||
adata: An AnnData object. | ||
layer: Key for the adata.layers dict in which the processed data will be saved. | ||
Default value is ``colocml_preprocessing``. | ||
If None, ``adata.X`` will be overwritten with the processed data. | ||
median_filter_size: 2-dimensional filter size for the median filtering that | ||
will be performed per ion image. | ||
quantile_threshold: Float between 0 and 1. The function computes the quantile value per | ||
ion image and all pixels below the quantile threshold will be set to 0. | ||
|
||
Returns: | ||
None. The processed data is saved in ``layer``. If layer is net to None, ``adata.X`` will | ||
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more.
|
||
be overwritten | ||
|
||
Raises: | ||
ValueError: If no annotations are available in ``adata.X``. | ||
|
||
""" | ||
|
||
# Extract image array from anndata: | ||
imarray = anndata_to_image_array(adata=adata) | ||
|
||
if len(imarray) == 0: | ||
raise ValueError("AnnData contains no annotations. ColocML preprocessing cannot be called.") | ||
|
||
# Median filtering | ||
imarray = ndimage.median_filter(imarray, size=(1, median_filter_size[0], median_filter_size[1])) | ||
|
||
# reshaping | ||
imarray = imarray.reshape((imarray.shape[0], -1)) | ||
|
||
# Quantile thresholding | ||
mask = imarray < np.percentile(imarray, q=quantile_threshold * 100, axis=1)[:, np.newaxis] | ||
imarray[mask] = 0 | ||
|
||
# Inserting of processed data into adata object | ||
if layer is None: | ||
adata.X = imarray.transpose() | ||
else: | ||
adata.layers[layer] = imarray.transpose() | ||
|
||
|
||
def colocalization(adata: AnnData, layer: Optional[str] = "coloc_ml_preprocessing"): | ||
""" | ||
Colocalization of ion images using the cosine similarity metric. | ||
|
||
In combination with the ``colocML_preprocessing`` function, this metric performed best in the | ||
colocML publication (https://doi.org/10.1093/bioinformatics/btaa085). | ||
|
||
It is recommended to call the the ``coloc_ml_preprocessing`` function beforehand. | ||
|
||
Args: | ||
adata: An AnnData object. | ||
layer: Key for ``adata.layer`` from which the ionimage_data for preprocessing taken. | ||
If ``None``, ``adata.X`` is used. ``coloc_ml_preprocessing`` will save the preprocessed | ||
data per default in ``adata.layer['coloc_ml_preprocessing']``. | ||
|
||
Returns: | ||
None. The processed data is saved in ``adata.varp['colocalization']``. | ||
|
||
Raises: | ||
ValueError: If layer is not found in adata.layers. | ||
""" | ||
|
||
# Select data | ||
if layer is None: | ||
data = np.array(adata.X).transpose() | ||
elif layer in adata.layers.keys(): | ||
data = np.array(adata.layers[layer]).transpose() | ||
else: | ||
raise ValueError(f"Layer `{layer}` not found in adata.layers.") | ||
|
||
# Compute colocalization | ||
coloc = _pairwise_cosine_similarity(data) | ||
|
||
# Save colocalization | ||
adata.varp[COLOCALIZATION] = coloc | ||
|
||
|
||
def _pairwise_cosine_similarity(data: np.ndarray) -> np.ndarray: | ||
|
||
# Divide image vectors by euclidean norm | ||
norm_data = data / np.linalg.norm(data, axis=1, keepdims=True) | ||
|
||
# Compute pairwise cosine similarity | ||
similarity_matrix = np.dot(norm_data, norm_data.T) | ||
|
||
return similarity_matrix |
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
|
@@ -28,3 +28,5 @@ class COL: | |
XY: Final = (X, Y) | ||
YX: Final = (Y, X) | ||
YXC: Final = (Y, X, C) | ||
|
||
COLOCALIZATION = "colocalization" |
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -0,0 +1,118 @@ | ||
from typing import TYPE_CHECKING, Tuple | ||
|
||
import numpy as np | ||
import pandas as pd | ||
import pytest | ||
from anndata import AnnData | ||
|
||
if TYPE_CHECKING: | ||
import _pytest.fixtures | ||
|
||
from metaspace_converter import anndata_to_image_array | ||
from metaspace_converter.colocalization import colocalization, coloc_ml_preprocessing | ||
from metaspace_converter.constants import COL, COLOCALIZATION, METASPACE_KEY, X, Y | ||
from metaspace_converter.to_anndata import all_image_pixel_coordinates | ||
|
||
|
||
@pytest.fixture | ||
def adata_dummy(request: "_pytest.fixtures.SubRequest") -> Tuple[AnnData, float]: | ||
num_ions = request.param.get("num_ions", 0) | ||
height = request.param.get("height", 0) | ||
width = request.param.get("width", 0) | ||
quantile_threshold = request.param.get("quantile_threshold", 0) | ||
shape = (num_ions, height, width) | ||
# Create a stack of ion images where every pixel has a different value | ||
ion_images_stack = np.arange(np.prod(shape)).reshape(shape) | ||
coordinates = all_image_pixel_coordinates(shape[1:]) | ||
# Create an AnnData with pixel coordinates | ||
obs = pd.DataFrame( | ||
{COL.ion_image_pixel_y: coordinates[:, 0], COL.ion_image_pixel_x: coordinates[:, 1]} | ||
) | ||
adata = AnnData( | ||
X=ion_images_stack.reshape((height * width, num_ions)), | ||
obs=obs, | ||
uns={METASPACE_KEY: {"image_size": {X: width, Y: height}}}, | ||
) | ||
return adata, quantile_threshold | ||
|
||
|
||
@pytest.mark.parametrize( | ||
"adata_dummy", | ||
[ | ||
dict(num_ions=5, height=10, width=10, quantile_threshold=0), | ||
dict(num_ions=0, height=10, width=10, quantile_threshold=0), | ||
dict(num_ions=5, height=3, width=4, quantile_threshold=0), | ||
dict(num_ions=5, height=3, width=4, quantile_threshold=0.5) | ||
], | ||
indirect=["adata_dummy"], | ||
) | ||
def test_colocml_preprocessing(adata_dummy): | ||
|
||
COLOCML_LAYER = "coloc_ml_preprocessing" | ||
|
||
adata, quantile_threshold = adata_dummy | ||
|
||
if adata.X.shape[1] == 0: | ||
with pytest.raises(ValueError) as e_info: | ||
coloc_ml_preprocessing( | ||
adata, median_filter_size=(3, 3), quantile_threshold=quantile_threshold, | ||
layer=COLOCML_LAYER | ||
) | ||
else: | ||
|
||
# Quantile thresholding and median filtering is tested separately | ||
if quantile_threshold == 0: | ||
|
||
# Test median filtering | ||
coloc_ml_preprocessing( | ||
adata, median_filter_size=(3, 3), quantile_threshold=quantile_threshold, | ||
layer=COLOCML_LAYER | ||
) | ||
|
||
actual = anndata_to_image_array(adata) | ||
|
||
# Layer exists | ||
assert COLOCML_LAYER in adata.layers.keys() | ||
|
||
# Layer sizes match | ||
assert adata.X.shape == adata.layers[COLOCML_LAYER].shape | ||
|
||
# median filter | ||
# Probe a single pixel of the resulting array and compute the expected median filter. | ||
x = 1 | ||
y = 1 | ||
width = adata.uns[METASPACE_KEY]["image_size"][X] | ||
print(actual[:, y - 1:y + 2, x - 1:x + 2]) | ||
expected_pixel_features = np.median(actual[:, y - 1:y + 2, x - 1:x + 2], axis=[1, 2]) | ||
actual_pixel_features = adata.layers[COLOCML_LAYER][y * width + x, :] | ||
np.testing.assert_allclose(actual_pixel_features, expected_pixel_features) | ||
|
||
# Quantile thresholding | ||
else: | ||
coloc_ml_preprocessing( | ||
adata, median_filter_size=(3, 3), quantile_threshold=quantile_threshold, | ||
layer=COLOCML_LAYER | ||
) | ||
quantile_value = quantile_threshold * adata.X.shape[0] | ||
actual_thresholded = np.sum(adata.layers[COLOCML_LAYER] == 0, axis=0) | ||
assert all(actual_thresholded <= quantile_value) | ||
|
||
|
||
|
||
|
||
|
||
def test_colocalization(): | ||
|
||
arr = np.array([[0.0, 1.0, 0.0], [0.0, 2.0, 0.0], [1.0, 0.0, 0.0]]).transpose() | ||
|
||
expected = np.array([[1.0, 1.0, 0.0], [1.0, 1.0, 0.0], [0.0, 0.0, 1.0]]) | ||
|
||
adata = AnnData(X=arr) | ||
|
||
colocalization(adata, layer=None) | ||
|
||
assert COLOCALIZATION in adata.varp.keys() | ||
|
||
coloc = adata.varp[COLOCALIZATION] | ||
|
||
assert np.all(coloc == expected) |
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
|
@@ -4,7 +4,7 @@ build-backend = "setuptools.build_meta" | |
|
||
[project] | ||
name = "metaspace-converter" | ||
version = "1.0.1" | ||
version = "1.1.0" | ||
authors = [ | ||
{name = "Tim Daniel Rose", email = "[email protected]"}, | ||
{name = "Andreas Eisenbarth", email = "[email protected]"}, | ||
|
Oops, something went wrong.
Add this suggestion to a batch that can be applied as a single commit.
This suggestion is invalid because no changes were made to the code.
Suggestions cannot be applied while the pull request is closed.
Suggestions cannot be applied while viewing a subset of changes.
Only one suggestion per line can be applied in a batch.
Add this suggestion to a batch that can be applied as a single commit.
Applying suggestions on deleted lines is not supported.
You must change the existing code in this line in order to create a valid suggestion.
Outdated suggestions cannot be applied.
This suggestion has been applied or marked resolved.
Suggestions cannot be applied from pending reviews.
Suggestions cannot be applied on multi-line comments.
Suggestions cannot be applied while the pull request is queued to merge.
Suggestion cannot be applied right now. Please check back later.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Also added
layer
argument to theanndata_to_image_array
function. With this also the processed images from layers can be extracted.