mirror of
https://github.com/chanzuckerberg/cellxgene.git
synced 2026-09-23 15:08:12 +08:00
Move cxgtool into CLI and modularize conversion functions (#1701)
This commit is contained in:
@@ -0,0 +1,131 @@
|
||||
from os import path
|
||||
|
||||
import click
|
||||
|
||||
from server.converters.h5ad_data_file import H5ADDataFile
|
||||
|
||||
|
||||
@click.command(
|
||||
name="convert",
|
||||
short_help="Converts an H5AD dataset to the CXG format.",
|
||||
help="Converts an H5AD dataset to the CXG format. The CXG format is a cellxgene-private data format "
|
||||
"that has performance and access characteristics amenable to a multi-dataset, multi-user serving "
|
||||
"environment. You will be able to launch the cellxgene using the `cellxgene launch` command as "
|
||||
"usually with the generated CXG file.",
|
||||
)
|
||||
@click.argument(
|
||||
"input-file",
|
||||
nargs=1,
|
||||
help="Path to the H5AD input file to be converted.",
|
||||
type=click.Path(exists=True, dir_okay=False),
|
||||
)
|
||||
@click.option(
|
||||
"-o",
|
||||
"--output-dir",
|
||||
help="Name of the output CXG directory. If not provided, will default to be the input filename with a "
|
||||
"CXG extension.",
|
||||
)
|
||||
@click.option(
|
||||
"-b",
|
||||
"--backed",
|
||||
help="When true, loads the H5AD in file backed mode. This will cause the conversion to be slower, "
|
||||
"but will use less memory.",
|
||||
default=False,
|
||||
show_default=True,
|
||||
is_flag=True,
|
||||
)
|
||||
@click.option(
|
||||
"-t",
|
||||
"--title",
|
||||
help="Human readable dataset title that will be included as metadata about the CXG file. If omitted, "
|
||||
"the dataset title will be the filename.",
|
||||
)
|
||||
@click.option(
|
||||
"-a",
|
||||
"--about",
|
||||
help="A fully qualified URL that provides more information about the dataset and will be included as "
|
||||
"metadata about the CXG file.",
|
||||
)
|
||||
@click.option(
|
||||
"-s",
|
||||
"--sparse-threshold",
|
||||
help="If the dataset's percent of non-zero values falls belows the specified threshold, then the X "
|
||||
"array of the dataset will be sparse. Since the default value is 0.0, the default will be to "
|
||||
"convert to dense array.",
|
||||
default=0.0,
|
||||
show_default=True,
|
||||
)
|
||||
@click.option("--obs-names",
|
||||
help="Name to a column in the obs dataframe that will be used as the index for the dataframe instead of "
|
||||
"the one designated by the dataframe generated-index.")
|
||||
@click.option("--var-names",
|
||||
help="Name to a column in the var dataframe that will be used as the index for the dataframe instead of "
|
||||
"the one designated by the dataframe generated-index.")
|
||||
@click.option(
|
||||
"--disable-custom-colors",
|
||||
help="When set, conversion process will not extract scanpy-compatible category colors from the H5AD file.",
|
||||
default=False,
|
||||
show_default=True,
|
||||
is_flag=True,
|
||||
)
|
||||
@click.option(
|
||||
"--disable-corpora-schema",
|
||||
"When set, conversion process will neither extract nor store Corpora schema information. See "
|
||||
"https://github.com/chanzuckerberg/corpora-data-portal/blob/main/backend/schema/corpora_schema.md for "
|
||||
"more information.",
|
||||
default=False,
|
||||
show_default=True,
|
||||
is_flag=True,
|
||||
)
|
||||
@click.option(
|
||||
"--overwrite",
|
||||
help="When set to true, will overwrite the output file if the output file already exists.",
|
||||
default=False,
|
||||
show_default=True,
|
||||
is_flag=True,
|
||||
)
|
||||
@click.option("-v", "--verbose", count=True)
|
||||
@click.help_option("--help", "-h", help="Show this message and exit.")
|
||||
def convert_to_cxg(
|
||||
input_file,
|
||||
output_directory,
|
||||
backed,
|
||||
title,
|
||||
about,
|
||||
sparse_threshold,
|
||||
obs_names,
|
||||
var_names,
|
||||
disable_custom_colors,
|
||||
disable_corpora_schema,
|
||||
should_overwrite,
|
||||
):
|
||||
"""
|
||||
Convert a dataset file into CXG.
|
||||
"""
|
||||
|
||||
h5ad_data_file = H5ADDataFile(input_file, backed, title, about, obs_names, var_names,
|
||||
use_corpora_schema=not disable_corpora_schema)
|
||||
|
||||
# Get the directory that will hold all the CXG files
|
||||
cxg_output_container = get_output_directory(input_file, output_directory, should_overwrite)
|
||||
|
||||
h5ad_data_file.to_cxg(cxg_output_container, sparse_threshold,
|
||||
convert_anndata_colors_to_cxg_colors=not disable_custom_colors)
|
||||
|
||||
|
||||
def get_output_directory(input_filename, output_directory, should_overwrite):
|
||||
"""
|
||||
Get the name of the CXG output directory to be created/populated during the dataset conversion.
|
||||
"""
|
||||
|
||||
if not path.isdir(output_directory) or (path.isdir(output_directory) and should_overwrite):
|
||||
if output_directory.endswith(".cxg"):
|
||||
return output_directory
|
||||
return output_directory + ".cxg"
|
||||
if path.isdir(output_directory) and not should_overwrite:
|
||||
raise click.BadParameter(
|
||||
f"Output directory {output_directory} already exists. If you'd like to overwrite, then run the command "
|
||||
f"with the --overwrite flag."
|
||||
)
|
||||
|
||||
return path.splitext(input_filename)[1] + ".cxg"
|
||||
@@ -8,7 +8,9 @@ import tiledb
|
||||
from flask import current_app
|
||||
|
||||
from server.common.annotations.annotations import Annotations
|
||||
from server.converters.cxgtool import sanitize_keys, generate_schema_hints_and_convert_value_types, cxg_dtype
|
||||
from server.common.errors import AnnotationCategoryNameError
|
||||
from server.common.utils.sanitization_utils import sanitize_values_in_list
|
||||
from server.common.utils.type_conversion_utils import get_dtypes_and_schemas_of_dataframe, get_dtype_of_array
|
||||
from server.db.cellxgene_orm import CellxGeneDataset, Annotation
|
||||
|
||||
|
||||
@@ -21,7 +23,12 @@ class AnnotationsHostedTileDB(Annotations):
|
||||
self.directory_path = directory_path
|
||||
|
||||
def check_category_names(self, df):
|
||||
sanitize_keys(df.keys().to_list(), False)
|
||||
original_category_names = df.keys().to_list()
|
||||
sanitized_category_names = set(sanitize_values_in_list(original_category_names).values())
|
||||
unsanitary_original_category_names = set(original_category_names).difference(sanitized_category_names)
|
||||
if unsanitary_original_category_names:
|
||||
raise AnnotationCategoryNameError(
|
||||
f"{unsanitary_original_category_names} are not valid category names, please resubmit")
|
||||
|
||||
def is_safe_collection_name(self, name):
|
||||
"""
|
||||
@@ -94,19 +101,19 @@ class AnnotationsHostedTileDB(Annotations):
|
||||
pass
|
||||
else:
|
||||
os.makedirs(uri, exist_ok=True)
|
||||
schema_hints, values = generate_schema_hints_and_convert_value_types(df)
|
||||
_, dataframe_schema_type_hints = get_dtypes_and_schemas_of_dataframe(df)
|
||||
|
||||
annotation = Annotation(
|
||||
tiledb_uri=uri,
|
||||
user_id=user_id,
|
||||
dataset_id=str(dataset_id),
|
||||
schema_hints=json.dumps(schema_hints)
|
||||
schema_hints=json.dumps(dataframe_schema_type_hints)
|
||||
)
|
||||
if not df.empty:
|
||||
self.check_category_names(df)
|
||||
# convert to tiledb datatypes
|
||||
for col in df:
|
||||
df[col] = df[col].astype(cxg_dtype(df[col]))
|
||||
df[col] = df[col].astype(get_dtype_of_array(df[col]))
|
||||
tiledb.from_pandas(uri, df)
|
||||
|
||||
self.db.session.add(annotation)
|
||||
|
||||
@@ -768,7 +768,7 @@ class DatasetConfig(BaseConfig):
|
||||
self.user_annotations__ontology__obo_location = dc["user_annotations"]["ontology"]["obo_location"]
|
||||
self.user_annotations__hosted_tiledb_array__db_uri = dc["user_annotations"]["hosted_tiledb_array"]["db_uri"]
|
||||
self.user_annotations__hosted_tiledb_array__hosted_file_directory = \
|
||||
dc["user_annotations"]["hosted_tiledb_array"]["hosted_file_directory"] # noqa E501
|
||||
dc["user_annotations"]["hosted_tiledb_array"]["hosted_file_directory"] # noqa E501
|
||||
|
||||
self.embeddings__names = dc["embeddings"]["names"]
|
||||
self.embeddings__enable_reembedding = dc["embeddings"]["enable_reembedding"]
|
||||
|
||||
@@ -9,6 +9,7 @@ import collections
|
||||
import json
|
||||
|
||||
from server.cli.upgrade import validate_version_str
|
||||
from server.common.utils.corpora_constants import CorporaConstants
|
||||
|
||||
|
||||
def corpora_get_versions_from_anndata(adata):
|
||||
@@ -56,26 +57,13 @@ def corpora_get_props_from_anndata(adata):
|
||||
if not version_is_supported:
|
||||
raise ValueError("Unsupported Corpora schema version")
|
||||
|
||||
required_simple_fields = [
|
||||
"version",
|
||||
"title",
|
||||
"layer_descriptions",
|
||||
"organism",
|
||||
"organism_ontology_term_id",
|
||||
"project_name",
|
||||
"project_description",
|
||||
]
|
||||
# Spec says some values encoded as JSON due to the inability of AnnData to store complex types.
|
||||
required_json_fields = ["contributors", "project_links"]
|
||||
optional_simple_fields = ["preprint_doi", "publication_doi", "default_embedding", "default_field", "tags"]
|
||||
|
||||
corpora_props = {}
|
||||
for key in required_simple_fields:
|
||||
for key in CorporaConstants.REQUIRED_SIMPLE_METADATA_FIELDS:
|
||||
if key not in adata.uns:
|
||||
raise KeyError(f"missing Corpora schema field {key}")
|
||||
corpora_props[key] = adata.uns[key]
|
||||
|
||||
for key in required_json_fields:
|
||||
for key in CorporaConstants.REQUIRED_JSON_ENCODED_METADATA_FIELD:
|
||||
if key not in adata.uns:
|
||||
raise KeyError(f"missing Corpora schema field {key}")
|
||||
try:
|
||||
@@ -83,7 +71,7 @@ def corpora_get_props_from_anndata(adata):
|
||||
except json.JSONDecodeError:
|
||||
raise json.JSONDecodeError(f"Corpora schema field {key} is expected to be a valid JSON string")
|
||||
|
||||
for key in optional_simple_fields:
|
||||
for key in CorporaConstants.OPTIONAL_SIMPLE_METADATA_FIELDS:
|
||||
if key in adata.uns:
|
||||
corpora_props[key] = adata.uns[key]
|
||||
|
||||
|
||||
@@ -0,0 +1,16 @@
|
||||
class CorporaConstants(object):
|
||||
REQUIRED_SIMPLE_METADATA_FIELDS = [
|
||||
"version",
|
||||
"title",
|
||||
"layer_descriptions",
|
||||
"organism",
|
||||
"organism_ontology_term_id",
|
||||
"project_name",
|
||||
"project_description",
|
||||
]
|
||||
|
||||
# The Corpora specification requires some values encoded as JSON due to the inability of AnnData to store complex
|
||||
# types.
|
||||
REQUIRED_JSON_ENCODED_METADATA_FIELD = ["contributors", "project_links"]
|
||||
|
||||
OPTIONAL_SIMPLE_METADATA_FIELDS = ["preprint_doi", "publication_doi", "default_embedding", "default_field", "tags"]
|
||||
@@ -0,0 +1,4 @@
|
||||
class CxgConstants(object):
|
||||
# The CXG container version number. Must be a semver string (major.minor.patch)
|
||||
# DO NOT UPDATE THIS WITHOUT ALSO UPDATING CXG SPECIFICATION.
|
||||
CXG_VERSION = "0.2.0"
|
||||
@@ -0,0 +1,179 @@
|
||||
import json
|
||||
|
||||
import numpy as np
|
||||
import tiledb
|
||||
|
||||
from server.common.utils.type_conversion_utils import get_dtype_of_array, get_dtype_and_schema_of_array
|
||||
|
||||
|
||||
def convert_dictionary_to_cxg_group(cxg_container, metadata_dict, group_metadata_name="cxg_group_metadata"):
|
||||
"""
|
||||
Saves the contents of the dictionary to the CXG output directory specified.
|
||||
|
||||
This function is primarily used to save metadata about a dataset to the CXG directory. At some point, tiledb will
|
||||
have support for metadata on groups at which point the utility of this function should be revisited. Until such
|
||||
feature exists, this function create an empty array and annotate that array.
|
||||
|
||||
For more information, visit https://github.com/TileDB-Inc/TileDB-Py/issues/254.
|
||||
"""
|
||||
|
||||
array_name = f"{cxg_container}/{group_metadata_name}"
|
||||
|
||||
# Because TileDB does not allow one to attach metadata directly to a CXG group, we need to have a workaround
|
||||
# where we create an empty array and attached the metadata onto to this empty array. Below we construct this empty
|
||||
# array.
|
||||
tiledb.from_numpy(array_name, np.zeros((1,)))
|
||||
|
||||
with tiledb.DenseArray(array_name, mode="w") as metadata_array:
|
||||
for key, value in metadata_dict.items():
|
||||
metadata_array.meta[key] = value
|
||||
|
||||
|
||||
def convert_dataframe_to_cxg_array(cxg_container, dataframe_name, dataframe, index_column_name, ctx):
|
||||
"""
|
||||
Saves the contents of the dataframe to the CXG output directory specified.
|
||||
|
||||
Current access patterns are oriented toward reading very large slices of the dataframe, one attribute at a time.
|
||||
Attribute data also tends to be (often) repetitive (bools, categories, strings). Given this, we use a large tile
|
||||
size (1000) and very aggressive compression levels.
|
||||
"""
|
||||
|
||||
def create_dataframe_array(array_name, dataframe):
|
||||
tiledb_filter = tiledb.FilterList(
|
||||
[
|
||||
# Attempt aggressive compression as many of these dataframes are very repetitive strings, bools and
|
||||
# other non-float data.
|
||||
tiledb.ZstdFilter(level=22),
|
||||
]
|
||||
)
|
||||
attrs = [
|
||||
tiledb.Attr(name=column, dtype=get_dtype_of_array(dataframe[column]), filters=tiledb_filter)
|
||||
for column in dataframe
|
||||
]
|
||||
domain = tiledb.Domain(
|
||||
tiledb.Dim(domain=(0, dataframe.shape[0] - 1), tile=min(dataframe.shape[0], 1000), dtype=np.uint32)
|
||||
)
|
||||
schema = tiledb.ArraySchema(
|
||||
domain=domain, sparse=False, attrs=attrs, cell_order="row-major", tile_order="row-major"
|
||||
)
|
||||
tiledb.DenseArray.create(array_name, schema)
|
||||
|
||||
array_name = f"{cxg_container}/{dataframe_name}"
|
||||
|
||||
create_dataframe_array(array_name, dataframe)
|
||||
|
||||
with tiledb.DenseArray(array_name, mode="w", ctx=ctx) as array:
|
||||
value = {}
|
||||
schema_hints = {}
|
||||
for column_name, column_values in dataframe.items():
|
||||
dtype, hints = get_dtype_and_schema_of_array(column_values)
|
||||
|
||||
value[column_name] = column_values.to_numpy(dtype=dtype)
|
||||
if hints:
|
||||
schema_hints.update({column_name: hints})
|
||||
|
||||
schema_hints.update({"index": index_column_name})
|
||||
array[:] = value
|
||||
array.meta["cxg_schema"] = json.dumps(schema_hints)
|
||||
|
||||
tiledb.consolidate(array_name, ctx=ctx)
|
||||
|
||||
|
||||
def convert_ndarray_to_cxg_dense_array(ndarray_name, ndarray, ctx):
|
||||
"""
|
||||
Saves contents of ndarray to the CXG output directory specified.
|
||||
|
||||
Generally this function is used to convert dataset embeddings. Because embeddings are typically accessed with
|
||||
very large slices (or all of the embedding), they do not benefit from overly aggressive compression due to their
|
||||
format. Given this, we use a large tile size (1000) but only default compression level.
|
||||
"""
|
||||
|
||||
def create_ndarray_array(ndarray_name, ndarray):
|
||||
filters = tiledb.FilterList([tiledb.ZstdFilter()])
|
||||
attrs = [tiledb.Attr(dtype=ndarray.dtype, filters=filters)]
|
||||
dimensions = [
|
||||
tiledb.Dim(
|
||||
domain=(0, ndarray.shape[dimension] - 1), tile=min(ndarray.shape[dimension], 1000), dtype=np.uint32
|
||||
)
|
||||
for dimension in range(ndarray.ndim)
|
||||
]
|
||||
domain = tiledb.Domain(*dimensions)
|
||||
schema = tiledb.ArraySchema(
|
||||
domain=domain, sparse=False, attrs=attrs, capacity=1_000_000, cell_order="row-major", tile_order="row-major"
|
||||
)
|
||||
tiledb.DenseArray.create(ndarray_name, schema)
|
||||
|
||||
create_ndarray_array(ndarray_name, ndarray)
|
||||
|
||||
with tiledb.DenseArray(ndarray_name, mode="w", ctx=ctx) as array:
|
||||
array[:] = ndarray
|
||||
|
||||
tiledb.consolidate(ndarray_name, ctx=ctx)
|
||||
|
||||
|
||||
def convert_matrix_to_cxg_array(
|
||||
matrix_name, matrix, encode_as_sparse_array, ctx, column_shift_for_sparse_encoding=None
|
||||
):
|
||||
"""
|
||||
Converts a numpy array matrix into a TileDB SparseArray of DenseArray based on whether `encode_as_sparse_array`
|
||||
is true or not. Note that when the matrix is encoded as a SparseArray, it only writes the values that are
|
||||
nonzero. This means that if you count the number of elements in the SparseArray, it will not equal the total
|
||||
number of elements in the matrix, only the number of nonzero elements.
|
||||
|
||||
Furthermore, if the `column_shift_for_sparse_encoding` matrix is not None, this function will subtract the sparse
|
||||
encoding from the original given matrix and as previously stated, only write the nonzero values to the TileDB
|
||||
SparseArray.
|
||||
"""
|
||||
|
||||
def create_matrix_array(matrix_name, number_of_rows, number_of_columns, encode_as_sparse_array):
|
||||
filters = tiledb.FilterList([tiledb.ZstdFilter()])
|
||||
attrs = [tiledb.Attr(dtype=np.float32, filters=filters)]
|
||||
if encode_as_sparse_array:
|
||||
domain = tiledb.Domain(
|
||||
tiledb.Dim(name="obs", domain=(0, number_of_rows - 1), tile=min(number_of_rows, 512), dtype=np.uint32),
|
||||
tiledb.Dim(
|
||||
name="var", domain=(0, number_of_columns - 1), tile=min(number_of_columns, 2048), dtype=np.uint32
|
||||
),
|
||||
)
|
||||
else:
|
||||
domain = tiledb.Domain(
|
||||
tiledb.Dim(name="obs", domain=(0, number_of_rows - 1), tile=min(number_of_rows, 50), dtype=np.uint32),
|
||||
tiledb.Dim(
|
||||
name="var", domain=(0, number_of_columns - 1), tile=min(number_of_columns, 100), dtype=np.uint32
|
||||
),
|
||||
)
|
||||
schema = tiledb.ArraySchema(
|
||||
domain=domain, sparse=encode_as_sparse_array, attrs=attrs, cell_order="row-major", tile_order="col-major"
|
||||
)
|
||||
if encode_as_sparse_array:
|
||||
tiledb.SparseArray.create(matrix_name, schema)
|
||||
else:
|
||||
tiledb.DenseArray.create(matrix_name, schema)
|
||||
|
||||
number_of_rows = matrix.shape[0]
|
||||
number_of_columns = matrix.shape[1]
|
||||
stride = min(int(np.power(10, np.around(np.log10(1e9 / number_of_columns)))), 10_000)
|
||||
|
||||
create_matrix_array(matrix_name, number_of_rows, number_of_columns, encode_as_sparse_array)
|
||||
|
||||
if encode_as_sparse_array:
|
||||
with tiledb.SparseArray(matrix_name, mode="w", ctx=ctx) as array:
|
||||
for start_row_index in range(0, number_of_rows, stride):
|
||||
end_row_index = min(start_row_index + stride, number_of_rows)
|
||||
matrix_subset = matrix[start_row_index:end_row_index, :]
|
||||
if not isinstance(matrix_subset, np.ndarray):
|
||||
matrix_subset = matrix_subset.toarray()
|
||||
if column_shift_for_sparse_encoding is not None:
|
||||
matrix_subset = matrix_subset - column_shift_for_sparse_encoding
|
||||
indices = np.nonzero(matrix_subset)
|
||||
trow = indices[0] + start_row_index
|
||||
array[trow, indices[1]] = matrix_subset[indices[0], indices[1]]
|
||||
|
||||
else:
|
||||
with tiledb.DenseArray(matrix_name, mode="w", ctx=ctx) as array:
|
||||
for start_row_index in range(0, number_of_rows, stride):
|
||||
end_row_index = min(start_row_index + stride, number_of_rows)
|
||||
matrix_subset = matrix[start_row_index:end_row_index, :]
|
||||
if not isinstance(matrix_subset, np.ndarray):
|
||||
matrix_subset = matrix_subset.toarray()
|
||||
array[start_row_index:end_row_index, :] = matrix_subset
|
||||
@@ -4,6 +4,17 @@ import numpy as np
|
||||
import pandas as pd
|
||||
|
||||
|
||||
def get_dtypes_and_schemas_of_dataframe(dataframe: pd.DataFrame):
|
||||
dtypes_by_column_name = {}
|
||||
schema_type_hints_by_column_name = {}
|
||||
|
||||
for column_name, column_values in dataframe.items():
|
||||
dtypes_by_column_name[column_name], schema_type_hints_by_column_name[column_name] = \
|
||||
get_dtype_and_schema_of_array(column_values)
|
||||
|
||||
return dtypes_by_column_name, schema_type_hints_by_column_name
|
||||
|
||||
|
||||
def get_dtype_of_array(array: pd.Series):
|
||||
return get_dtype_and_schema_of_array(array)[0]
|
||||
|
||||
|
||||
@@ -1,669 +0,0 @@
|
||||
"""
|
||||
This program converts an [AnnData H5AD](https://anndata.readthedocs.io/en/stable/)
|
||||
into a cellxgene TileDB structure, aka a [CXG](../../dev_docs/cxg.md).
|
||||
|
||||
IF YOU UPDATE THIS FILE, IN ANY WAY THAT MODIFIES THE CXG FORMAT or CONTENTS,
|
||||
YOU MUST UPDATE THE CXG SPECIFICATION and VERSION NUMBER.
|
||||
"""
|
||||
import re
|
||||
import anndata
|
||||
import tiledb
|
||||
import argparse
|
||||
import numpy as np
|
||||
from os.path import splitext, basename
|
||||
import json
|
||||
from scipy.stats import mode
|
||||
|
||||
from server.common.colors import convert_anndata_category_colors_to_cxg_category_colors
|
||||
from server.common.errors import ColorFormatException, AnnotationCategoryNameError
|
||||
from server.common.corpora import (
|
||||
corpora_get_props_from_anndata,
|
||||
corpora_get_versions_from_anndata,
|
||||
corpora_is_version_supported,
|
||||
)
|
||||
|
||||
|
||||
# the CXG container version number. Must be a semver string (major.minor.patch)
|
||||
# DO NOT UPDATE THIS WITHOUT ALSO UPDATING THE CXG SPECIFICATION.
|
||||
CXG_VERSION = "0.2.0"
|
||||
|
||||
# log_level must have a default
|
||||
log_level = 3
|
||||
|
||||
|
||||
def log(level, *args):
|
||||
global log_level
|
||||
if log_level and level <= log_level:
|
||||
print(*args)
|
||||
|
||||
|
||||
def main():
|
||||
parser = argparse.ArgumentParser()
|
||||
parser.add_argument("h5ad", nargs="?", help="H5AD file name")
|
||||
parser.add_argument(
|
||||
"--backed", action="store_true", help="loaded in file backed mode. Will be slower, but use less memory."
|
||||
)
|
||||
parser.add_argument(
|
||||
"--disable-custom-colors",
|
||||
action="store_true",
|
||||
default=False,
|
||||
help="Do not extract scanpy-compatible category colors from h5ad file.",
|
||||
)
|
||||
parser.add_argument(
|
||||
"--obs-names", help="Name of annotation to use for observations. If not specified, will use the obs index."
|
||||
)
|
||||
parser.add_argument(
|
||||
"--var-names", help="Name of annotation to use for variables. If not specified, will use the var index."
|
||||
)
|
||||
parser.add_argument("--verbose", "-v", action="count", default=0, help="verbose output")
|
||||
parser.add_argument("--title", help="Human readable dataset title. If omitted, will use filename")
|
||||
parser.add_argument(
|
||||
"--about",
|
||||
metavar="<URL>",
|
||||
help="URL providing more information about the dataset (hint: must be a fully specified absolute URL).",
|
||||
)
|
||||
parser.add_argument("--out", "--output", "-o", help="output CXG file name")
|
||||
parser.add_argument(
|
||||
"--sparse-threshold",
|
||||
"-s",
|
||||
type=float,
|
||||
default=0.0, # force dense by default
|
||||
help="The X array will be sparse if the percent of non-zeros falls below this value",
|
||||
)
|
||||
parser.add_argument(
|
||||
"--disable-corpora",
|
||||
action="store_true",
|
||||
default=False,
|
||||
help="Disable extraction and storing of Corpora schema information.",
|
||||
)
|
||||
args = parser.parse_args()
|
||||
|
||||
global log_level
|
||||
log_level = args.verbose
|
||||
|
||||
adata = anndata.read_h5ad(args.h5ad, backed="r" if args.backed else None)
|
||||
log(1, f"{basename(args.h5ad)} loaded...")
|
||||
|
||||
basefname = splitext(basename(args.h5ad))[0]
|
||||
out = args.out if args.out is not None else basefname
|
||||
container = out if splitext(out)[1] == ".cxg" else out + ".cxg"
|
||||
|
||||
corpora_props = load_corpora_props(args, adata) if not args.disable_corpora else None
|
||||
cxg_group_metadata = create_cxg_group_metadata(
|
||||
adata,
|
||||
basefname,
|
||||
title=args.title,
|
||||
about=args.about,
|
||||
corpora_props=corpora_props,
|
||||
extract_colors=not args.disable_custom_colors,
|
||||
)
|
||||
|
||||
write_cxg(
|
||||
adata,
|
||||
container,
|
||||
cxg_group_metadata=cxg_group_metadata,
|
||||
var_names=args.var_names,
|
||||
obs_names=args.obs_names,
|
||||
sparse_threshold=args.sparse_threshold,
|
||||
)
|
||||
|
||||
log(1, "done")
|
||||
|
||||
|
||||
def write_cxg(adata, container, cxg_group_metadata, var_names=None, obs_names=None, sparse_threshold=5.0):
|
||||
if not adata.var.index.is_unique:
|
||||
raise ValueError("Variable index is not unique - unable to convert.")
|
||||
if not adata.obs.index.is_unique:
|
||||
raise ValueError("Observation index is not unique - unable to convert.")
|
||||
|
||||
"""
|
||||
TileDB bug TileDB-Inc/TileDB#1575 requires that we sanitize all column names
|
||||
prior to saving. This can be reverted when the bug is fixed.
|
||||
"""
|
||||
log(0, "Warning: sanitizing all dataframe column names.")
|
||||
clean_all_column_names(adata)
|
||||
|
||||
ctx = tiledb.Ctx(
|
||||
{
|
||||
"sm.num_reader_threads": 32,
|
||||
"sm.num_writer_threads": 32,
|
||||
"sm.consolidation.buffer_size": 1 * 1024 * 1024 * 1024,
|
||||
}
|
||||
)
|
||||
|
||||
tiledb.group_create(container, ctx=ctx)
|
||||
log(1, f"\t...group created, with name {container}")
|
||||
|
||||
# dataset metadata
|
||||
save_metadata(container, cxg_group_metadata)
|
||||
log(1, "\t...dataset metadata saved")
|
||||
|
||||
# var/gene dataframe
|
||||
save_dataframe(container, "var", adata.var, var_names, ctx=ctx)
|
||||
log(1, "\t...var dataframe created")
|
||||
|
||||
# obs/cell dataframe
|
||||
save_dataframe(container, "obs", adata.obs, obs_names, ctx=ctx)
|
||||
log(1, "\t...obs dataframe created")
|
||||
|
||||
# embeddings
|
||||
e_container = f"{container}/emb"
|
||||
tiledb.group_create(e_container, ctx=ctx)
|
||||
save_embeddings(e_container, adata, ctx)
|
||||
log(1, "\t...embeddings created")
|
||||
|
||||
# X matrix
|
||||
save_X(container, adata.X, ctx, sparse_threshold)
|
||||
log(1, "\t...X created")
|
||||
|
||||
|
||||
"""
|
||||
TODO: the code used to handle type inferencing should not be duplicated between
|
||||
this tool and the server/common/utils code. When this tool is merged into
|
||||
the cellxgene CLI, consolidate.
|
||||
"""
|
||||
|
||||
|
||||
def dtype_to_schema(dtype):
|
||||
if dtype == np.float32:
|
||||
return (np.float32, {})
|
||||
elif dtype == np.int32:
|
||||
return (np.int32, {})
|
||||
elif dtype == np.bool_:
|
||||
return (np.uint8, {"type": "boolean"})
|
||||
elif dtype == np.str:
|
||||
return (np.unicode, {"type": "string"})
|
||||
elif dtype == "category":
|
||||
typ, hint = cxg_type(dtype.categories)
|
||||
return (typ, {"type": "categorical", "categories": dtype.categories.tolist()})
|
||||
else:
|
||||
raise TypeError(f"Annotations of type {dtype} are unsupported.")
|
||||
|
||||
|
||||
def _can_cast_to_float32(array):
|
||||
if array.dtype.kind == "f":
|
||||
# force downcast for all floats
|
||||
return True
|
||||
return False
|
||||
|
||||
|
||||
def _can_cast_to_int32(array):
|
||||
if array.dtype.kind in ["i", "u"]:
|
||||
if np.can_cast(array.dtype, np.int32):
|
||||
return True
|
||||
ii32 = np.iinfo(np.int32)
|
||||
if array.min() >= ii32.min and array.max() <= ii32.max:
|
||||
return True
|
||||
return False
|
||||
|
||||
|
||||
def cxg_type(array):
|
||||
try:
|
||||
return dtype_to_schema(array.dtype)
|
||||
except TypeError:
|
||||
dtype = array.dtype
|
||||
data_kind = dtype.kind
|
||||
if _can_cast_to_float32(array):
|
||||
return (np.float32, {})
|
||||
elif _can_cast_to_int32(array):
|
||||
return (np.int32, {})
|
||||
elif data_kind == "O" and dtype == "object":
|
||||
return (np.unicode, {"type": "string"})
|
||||
else:
|
||||
raise TypeError(f"Annotations of type {dtype} are unsupported.")
|
||||
|
||||
|
||||
def cxg_dtype(array):
|
||||
return cxg_type(array)[0]
|
||||
|
||||
|
||||
def create_dataframe(name, df, ctx):
|
||||
"""
|
||||
Current access patterns are oriented toward reading very large slices of
|
||||
the dataframe, one attribute at a time. Attribute data also tends to be
|
||||
(often) repetitive (bools, categories, strings).
|
||||
Given this, we use:
|
||||
* a large tile size (1000)
|
||||
* very aggressive compression levels
|
||||
"""
|
||||
filter = tiledb.FilterList(
|
||||
[
|
||||
# attempt aggressive compression as many of these dataframes are very repetitive
|
||||
# strings, bools and other non-float data.
|
||||
tiledb.ZstdFilter(level=22),
|
||||
]
|
||||
)
|
||||
attrs = [tiledb.Attr(name=col, dtype=cxg_dtype(df[col]), filters=filter) for col in df]
|
||||
domain = tiledb.Domain(tiledb.Dim(domain=(0, df.shape[0] - 1), tile=min(df.shape[0], 1000), dtype=np.uint32))
|
||||
schema = tiledb.ArraySchema(
|
||||
domain=domain, sparse=False, attrs=attrs, cell_order="row-major", tile_order="row-major"
|
||||
)
|
||||
tiledb.DenseArray.create(name, schema)
|
||||
|
||||
|
||||
def create_unique_column_name(df_cols, col_name_prefix):
|
||||
"""
|
||||
given the columns of a dataframe, and a name prefix, return a column name which
|
||||
does not exist in the dataframe, AND which is prefixed by `prefix`
|
||||
|
||||
The approach is to append a numeric suffix, starting at zero and increasing by
|
||||
one, until an unused name is found (eg, prefix_0, prefix_1, ...).
|
||||
"""
|
||||
suffix = 0
|
||||
while f"{col_name_prefix}{suffix}" in df_cols:
|
||||
suffix += 1
|
||||
return f"{col_name_prefix}{suffix}"
|
||||
|
||||
|
||||
def alias_index_col(df, df_name, index_col_name):
|
||||
"""
|
||||
We rely in the existance of a unique, human-readable index for
|
||||
any dataframe (eg, var is typically gene name, obs the cell name).
|
||||
The user can specify these via the --obs-names and --var-names config.
|
||||
If they are not specified, use the existing index to create them, giving
|
||||
the resulting column a unique name (eg, "name").
|
||||
|
||||
In both cases, enforce that the result is unique, and communicate the
|
||||
index column name via the 'index' field in the schema hints.
|
||||
"""
|
||||
if index_col_name is None:
|
||||
if not df.index.is_unique:
|
||||
raise KeyError(
|
||||
f"Values in {df_name}.index must be unique. "
|
||||
"Please prepare data to contain unique index values, or specify an "
|
||||
"alternative with --{ax_name}-name."
|
||||
)
|
||||
index_col_name = create_unique_column_name(df.columns, "name_")
|
||||
# turn the index into a normal column
|
||||
df.rename_axis(index_col_name, inplace=True)
|
||||
df.reset_index(inplace=True)
|
||||
|
||||
elif index_col_name in df.columns:
|
||||
# User has specified alternative column for unique names, and it exists
|
||||
if not df[index_col_name].is_unique:
|
||||
raise KeyError(
|
||||
f"Values in {df_name}.{index_col_name} must be unique. Please prepare data to contain unique values."
|
||||
)
|
||||
|
||||
else:
|
||||
raise KeyError(f"Annotation {index_col_name}, specified in --{df_name}-name, does not exist.")
|
||||
|
||||
return (df, index_col_name)
|
||||
|
||||
|
||||
def generate_schema_hints_and_convert_value_types(df):
|
||||
value = {}
|
||||
schema_hints = {}
|
||||
for k, v in df.items():
|
||||
dtype, hints = cxg_type(v)
|
||||
value[k] = v.to_numpy(dtype=dtype)
|
||||
if hints:
|
||||
schema_hints.update({k: hints})
|
||||
return schema_hints, value
|
||||
|
||||
|
||||
def save_dataframe(container, name, df, index_col_name, ctx):
|
||||
A_name = f"{container}/{name}"
|
||||
(df, index_col_name) = alias_index_col(df, name, index_col_name)
|
||||
create_dataframe(A_name, df, ctx=ctx)
|
||||
with tiledb.DenseArray(A_name, mode="w", ctx=ctx) as A:
|
||||
schema_hints, value = generate_schema_hints_and_convert_value_types(df)
|
||||
schema_hints.update({"index": index_col_name})
|
||||
# convert all values in all cols to a numpy version of cxg datatypes,
|
||||
# then store the contents in the tiledb array A
|
||||
A[:] = value
|
||||
A.meta["cxg_schema"] = json.dumps(schema_hints)
|
||||
|
||||
tiledb.consolidate(A_name, ctx=ctx)
|
||||
|
||||
|
||||
def create_emb(e_name, emb):
|
||||
"""
|
||||
Embeddings are typically accessed with very large slices (or all of the embedding),
|
||||
and do not benefit from overly aggressive compression due to their format. Given
|
||||
this, we use:
|
||||
* large tile size (1000)
|
||||
* default compression level
|
||||
"""
|
||||
filters = tiledb.FilterList([tiledb.ZstdFilter()])
|
||||
attrs = [tiledb.Attr(dtype=emb.dtype, filters=filters)]
|
||||
dims = []
|
||||
for d in range(emb.ndim):
|
||||
shape = emb.shape
|
||||
dims.append(tiledb.Dim(domain=(0, shape[d] - 1), tile=min(shape[d], 1000), dtype=np.uint32))
|
||||
domain = tiledb.Domain(*dims)
|
||||
schema = tiledb.ArraySchema(
|
||||
domain=domain, sparse=False, attrs=attrs, capacity=1_000_000, cell_order="row-major", tile_order="row-major"
|
||||
)
|
||||
tiledb.DenseArray.create(e_name, schema)
|
||||
|
||||
|
||||
def is_valid_embedding(adata, name, arr):
|
||||
""" return True if this layout data is a valid array for front-end presentation:
|
||||
* ndarray, with shape (n_obs, >= 2), dtype float/int/uint
|
||||
* follows ScanPy embedding naming conventions
|
||||
* with all values finite or NaN (no +Inf or -Inf)
|
||||
"""
|
||||
is_valid = type(name) == str and name.startswith("X_") and len(name) > 2
|
||||
is_valid = is_valid and type(arr) == np.ndarray and arr.dtype.kind in "fiu"
|
||||
is_valid = is_valid and arr.shape[0] == adata.n_obs and arr.shape[1] >= 2
|
||||
is_valid = is_valid and not np.any(np.isinf(arr)) and not np.all(np.isnan(arr))
|
||||
return is_valid
|
||||
|
||||
|
||||
def save_embeddings(container, adata, ctx):
|
||||
for (name, value) in adata.obsm.items():
|
||||
if is_valid_embedding(adata, name, value):
|
||||
e_name = f"{container}/{name[2:]}"
|
||||
create_emb(e_name, value)
|
||||
with tiledb.DenseArray(e_name, mode="w", ctx=ctx) as A:
|
||||
A[:] = value
|
||||
tiledb.consolidate(e_name, ctx=ctx)
|
||||
log(1, f"\t\t...{name} embedding created")
|
||||
|
||||
|
||||
def create_X(X_name, shape, is_sparse):
|
||||
"""
|
||||
The X matrix is accessed in both row and column oriented patterns, depending on the
|
||||
particular operation. Because of the data type, default compression works best.
|
||||
The tile size, (50, 100) for dense, and (512,2048) for sparse,
|
||||
and global layout (row/col) was chosen empirically, by benchmarking
|
||||
the current cellxgene backend.
|
||||
"""
|
||||
filters = tiledb.FilterList([tiledb.ZstdFilter()])
|
||||
attrs = [tiledb.Attr(dtype=np.float32, filters=filters)]
|
||||
if is_sparse:
|
||||
domain = tiledb.Domain(
|
||||
tiledb.Dim(name="obs", domain=(0, shape[0] - 1), tile=min(shape[0], 512), dtype=np.uint32),
|
||||
tiledb.Dim(name="var", domain=(0, shape[1] - 1), tile=min(shape[1], 2048), dtype=np.uint32),
|
||||
)
|
||||
else:
|
||||
domain = tiledb.Domain(
|
||||
tiledb.Dim(name="obs", domain=(0, shape[0] - 1), tile=min(shape[0], 50), dtype=np.uint32),
|
||||
tiledb.Dim(name="var", domain=(0, shape[1] - 1), tile=min(shape[1], 100), dtype=np.uint32),
|
||||
)
|
||||
schema = tiledb.ArraySchema(
|
||||
domain=domain, sparse=is_sparse, attrs=attrs, cell_order="row-major", tile_order="col-major"
|
||||
)
|
||||
if is_sparse:
|
||||
tiledb.SparseArray.create(X_name, schema)
|
||||
else:
|
||||
tiledb.DenseArray.create(X_name, schema)
|
||||
|
||||
|
||||
def evaluate_for_sparse_encoding(xdata, sparse_threshold):
|
||||
"""
|
||||
This function determines if the X matrix has a sparsity below the sparse_threshold.
|
||||
This function also returns the number of non-zeros encountered and number
|
||||
of elements evaluated. This function may return before evaluating the whole X matrix
|
||||
if it can be determined that X is not sparse enough.
|
||||
"""
|
||||
shape = xdata.shape
|
||||
stride = min(int(np.power(10, np.around(np.log10(1e9 / shape[1])))), 10_000)
|
||||
nnz = 0
|
||||
maxnnz = int(shape[0] * shape[1] * sparse_threshold / 100)
|
||||
for row in range(0, shape[0], stride):
|
||||
lim = min(row + stride, shape[0])
|
||||
a = xdata[row:lim, :]
|
||||
if type(a) is not np.ndarray:
|
||||
a = a.toarray()
|
||||
nnz += np.count_nonzero(a)
|
||||
if nnz > maxnnz:
|
||||
return (False, nnz, lim * shape[1])
|
||||
log(2, "\t...rows", lim, "of", shape[0], "nnz", nnz, "nnz percent %5.2f%%" % (100 * nnz / (lim * shape[1])))
|
||||
|
||||
is_sparse = (100.0 * nnz / (shape[0] * shape[1])) < sparse_threshold
|
||||
return (is_sparse, nnz, shape[0] * shape[1])
|
||||
|
||||
|
||||
def evaluate_for_sparse_column_shift_encoding(xdata, sparse_threshold):
|
||||
"""Column shift encoding works by taking the most common value in each column, then
|
||||
subtracting that value from each element of the column. If each column mostly contains
|
||||
its most common value, then the resulting matrix can be very sparse.
|
||||
|
||||
This function determines if column shift encoding can be used to transform
|
||||
the X matrix into a sparse matrix with a sparsity below the sparse_threshold.
|
||||
If so, return the col_shift array that stores this encoding.
|
||||
This function also returns the number of non-zeros encountered and number
|
||||
of elements evaluated. This function may return before evaluating the whole X matrix
|
||||
if it can be determined that X cannot benefit from column shift encoding.
|
||||
"""
|
||||
shape = xdata.shape
|
||||
stride = max(1, 128_000_000 // shape[0])
|
||||
col_shift = np.zeros(shape[1])
|
||||
nnz = 0
|
||||
maxnnz = int(shape[0] * shape[1] * sparse_threshold / 100)
|
||||
for col in range(0, shape[1], stride):
|
||||
lim = min(col + stride, shape[1])
|
||||
a = xdata[:, col:lim]
|
||||
if type(a) is not np.ndarray:
|
||||
a = a.toarray()
|
||||
m = mode(a)
|
||||
col_shift[col:lim] = m.mode
|
||||
nnz += shape[0] * (lim - col) - np.sum(m.count)
|
||||
if nnz > maxnnz:
|
||||
return (None, nnz, shape[0] * lim)
|
||||
log(2, "\t...cols", lim, "of", shape[1], "nnz", nnz, "nnz percent %5.2f%%" % (100 * nnz / (lim * shape[0])))
|
||||
|
||||
is_sparse = (100.0 * nnz / (shape[0] * shape[1])) < sparse_threshold
|
||||
return (col_shift if is_sparse else None, nnz, shape[0] * shape[1])
|
||||
|
||||
|
||||
def save_X(container, xdata, ctx, sparse_threshold, expect_sparse=False):
|
||||
# Save X count matrix
|
||||
X_name = f"{container}/X"
|
||||
|
||||
shape = xdata.shape
|
||||
log(1, "\t...shape:", str(shape))
|
||||
|
||||
col_shift = None
|
||||
if sparse_threshold == 100:
|
||||
is_sparse = True
|
||||
elif sparse_threshold == 0:
|
||||
is_sparse = False
|
||||
else:
|
||||
is_sparse, nnz, nelem = evaluate_for_sparse_encoding(xdata, sparse_threshold)
|
||||
percent = 100.0 * nnz / nelem
|
||||
if nelem != shape[0] * shape[1]:
|
||||
log(1, "\t...sparse=", is_sparse, "non-zeros percent (estimate): %6.2f" % percent)
|
||||
else:
|
||||
log(1, "\t...sparse=", is_sparse, "non-zeros:", nnz, "percent: %6.2f" % percent)
|
||||
|
||||
is_sparse = percent < sparse_threshold
|
||||
if not is_sparse:
|
||||
col_shift, nnz, nelem = evaluate_for_sparse_column_shift_encoding(xdata, sparse_threshold)
|
||||
is_sparse = col_shift is not None
|
||||
percent = 100.0 * nnz / nelem
|
||||
if nelem != shape[0] * shape[1]:
|
||||
log(1, "\t...sparse=", is_sparse, "col shift non-zeros percent (estimate): %6.2f" % percent)
|
||||
else:
|
||||
log(1, "\t...sparse=", is_sparse, "col shift non-zeros:", nnz, "percent: %6.2f" % percent)
|
||||
|
||||
if expect_sparse is True and is_sparse is False:
|
||||
return False
|
||||
|
||||
create_X(X_name, shape, is_sparse)
|
||||
stride = min(int(np.power(10, np.around(np.log10(1e9 / shape[1])))), 10_000)
|
||||
if is_sparse:
|
||||
if col_shift is not None:
|
||||
log(1, "\t...output X as sparse matrix with column shift encoding")
|
||||
X_col_shift_name = f"{container}/X_col_shift"
|
||||
filters = tiledb.FilterList([tiledb.ZstdFilter()])
|
||||
attrs = [tiledb.Attr(dtype=np.float32, filters=filters)]
|
||||
domain = tiledb.Domain(tiledb.Dim(domain=(0, shape[1] - 1), tile=min(shape[1], 5000), dtype=np.uint32))
|
||||
schema = tiledb.ArraySchema(domain=domain, attrs=attrs)
|
||||
tiledb.DenseArray.create(X_col_shift_name, schema)
|
||||
with tiledb.DenseArray(X_col_shift_name, mode="w", ctx=ctx) as X_col_shift:
|
||||
X_col_shift[:] = col_shift
|
||||
tiledb.consolidate(X_col_shift_name, ctx=ctx)
|
||||
else:
|
||||
log(1, "\t...output X as sparse matrix")
|
||||
|
||||
with tiledb.SparseArray(X_name, mode="w", ctx=ctx) as X:
|
||||
nnz = 0
|
||||
for row in range(0, shape[0], stride):
|
||||
lim = min(row + stride, shape[0])
|
||||
a = xdata[row:lim, :]
|
||||
if type(a) is not np.ndarray:
|
||||
a = a.toarray()
|
||||
if col_shift is not None:
|
||||
a = a - col_shift
|
||||
indices = np.nonzero(a)
|
||||
trow = indices[0] + row
|
||||
nnz += indices[0].shape[0]
|
||||
X[trow, indices[1]] = a[indices[0], indices[1]]
|
||||
log(2, "\t...rows", lim, "of", shape[0], "nnz", nnz, "sparse", nnz / (lim * shape[1]))
|
||||
|
||||
else:
|
||||
log(1, "\t...output X as dense matrix")
|
||||
with tiledb.DenseArray(X_name, mode="w", ctx=ctx) as X:
|
||||
for row in range(0, shape[0], stride):
|
||||
lim = min(row + stride, shape[0])
|
||||
a = xdata[row:lim, :]
|
||||
if type(a) is not np.ndarray:
|
||||
a = a.toarray()
|
||||
X[row:lim, :] = a
|
||||
log(2, "\t...rows", row, "to", lim)
|
||||
|
||||
tiledb.consolidate(X_name, ctx=ctx)
|
||||
if hasattr(tiledb, "vacuum"):
|
||||
tiledb.vacuum(X_name)
|
||||
|
||||
return is_sparse
|
||||
|
||||
|
||||
def save_metadata(container, metadata_dict):
|
||||
"""
|
||||
Save all dataset-wide metadata. This includes:
|
||||
* CXG version
|
||||
* dataset metadata, such as title and about link.
|
||||
|
||||
Longer term, tiledb will have support for metadata on groups. Until
|
||||
such feature exists, create an empty array and annotate that array.
|
||||
|
||||
https://github.com/TileDB-Inc/TileDB-Py/issues/254
|
||||
"""
|
||||
a_name = f"{container}/cxg_group_metadata"
|
||||
with tiledb.from_numpy(a_name, np.zeros((1,))) as A:
|
||||
pass
|
||||
with tiledb.DenseArray(a_name, mode="w") as A:
|
||||
for k, v in metadata_dict.items():
|
||||
A.meta[k] = v
|
||||
|
||||
|
||||
def load_corpora_props(args, adata):
|
||||
versions = corpora_get_versions_from_anndata(adata)
|
||||
if versions is None:
|
||||
return None
|
||||
|
||||
[corpora_schema_version, corpora_encoding_version] = versions
|
||||
corpora_props = corpora_get_props_from_anndata(adata)
|
||||
version_is_supported = corpora_is_version_supported(corpora_schema_version, corpora_encoding_version)
|
||||
if not version_is_supported or not corpora_props:
|
||||
log(0, "ERROR: Unknown source file schema version is unsupported")
|
||||
raise ValueError("Unsupported Corpora schema version")
|
||||
|
||||
log(1, "FYI, file appears to be encoded using Corpora schema standards...")
|
||||
if args.title is not None or args.about is not None:
|
||||
log(0, "Warning: explicit specification of --title or --about will override Corpora schema fields.")
|
||||
|
||||
return corpora_props
|
||||
|
||||
|
||||
def create_cxg_group_metadata(adata, basefname, title=None, about=None, corpora_props=None, extract_colors=True):
|
||||
|
||||
if corpora_props is not None:
|
||||
# clobber encoding version to be OUR version, not the source H5AD encoding
|
||||
corpora_props["version"].update({"corpora_encoding_version": CXG_VERSION})
|
||||
corpora_project_links = corpora_props.get("project_links", [])
|
||||
corpora_about_link = next(
|
||||
(link for link in corpora_project_links if (link.get("link_type", None) == "SUMMARY")), {}
|
||||
)
|
||||
else:
|
||||
corpora_about_link = {}
|
||||
|
||||
title = title or corpora_about_link.get("link_name", basefname)
|
||||
about = about or corpora_about_link.get("link_url")
|
||||
|
||||
cxg_group_metadata = {"cxg_version": CXG_VERSION, "cxg_properties": json.dumps({"title": title, "about": about})}
|
||||
if corpora_props is not None:
|
||||
cxg_group_metadata.update({"corpora": json.dumps(corpora_props)})
|
||||
|
||||
if extract_colors:
|
||||
try:
|
||||
cxg_group_metadata["cxg_category_colors"] = json.dumps(
|
||||
convert_anndata_category_colors_to_cxg_category_colors(adata)
|
||||
)
|
||||
except ColorFormatException:
|
||||
log(
|
||||
0,
|
||||
"Warning: failed to extract colors from h5ad file! "
|
||||
"Fix the h5ad file or rerun with --disable-custom-colors. See help for details.",
|
||||
)
|
||||
|
||||
return cxg_group_metadata
|
||||
|
||||
|
||||
def sanitize_keys(keys, update_keys=True):
|
||||
"""
|
||||
We need names to be safe to use as attribute names in tiledb. See:
|
||||
TileDB-Inc/TileDB#1575
|
||||
TileDB-Inc/TileDB-Py#294
|
||||
This can be entirely removed once they add proper escaping.
|
||||
|
||||
Args: list of keys
|
||||
Returns: dict of {old_key: new_key, ...}
|
||||
|
||||
Returned new keys will be both safe and unique.
|
||||
|
||||
Masking out [~/.] and anything outside the ASCII range.
|
||||
"""
|
||||
p = re.compile(r"[^ -\.0-\[\]-\}]")
|
||||
clean_keys = {k: p.sub("_", k) for k in keys}
|
||||
|
||||
used_keys = set()
|
||||
clean_unique_keys = {}
|
||||
for k, v in clean_keys.items():
|
||||
if v not in used_keys:
|
||||
used_keys.add(v)
|
||||
clean_unique_keys[k] = v
|
||||
continue
|
||||
|
||||
# else, needs deduping.
|
||||
counter = 1
|
||||
while True:
|
||||
candidate_name = v + "-" + str(counter)
|
||||
if candidate_name not in used_keys:
|
||||
used_keys.add(candidate_name)
|
||||
clean_unique_keys[k] = candidate_name
|
||||
break
|
||||
counter += 1
|
||||
|
||||
for k, v, in clean_unique_keys.items():
|
||||
if k != v:
|
||||
if update_keys is False:
|
||||
raise AnnotationCategoryNameError(f"{k} not a valid category name, please resubmit")
|
||||
log(1, f"Renaming {k} to {v}")
|
||||
return clean_unique_keys
|
||||
|
||||
|
||||
def sanitize_df(df):
|
||||
df.rename(columns=sanitize_keys(df.keys().tolist()), inplace=True)
|
||||
|
||||
|
||||
def sanitize_mapping(mapping):
|
||||
clean_keys = sanitize_keys([k for k in mapping.keys()])
|
||||
for old_key, new_key in clean_keys.items():
|
||||
if old_key != new_key:
|
||||
mapping[new_key] = mapping[old_key]
|
||||
del mapping[old_key]
|
||||
|
||||
|
||||
def clean_all_column_names(adata):
|
||||
sanitize_df(adata.obs)
|
||||
sanitize_df(adata.var)
|
||||
sanitize_mapping(adata.obsm)
|
||||
|
||||
|
||||
if __name__ == "__main__":
|
||||
main()
|
||||
@@ -0,0 +1,250 @@
|
||||
import json
|
||||
import logging
|
||||
from os import path
|
||||
|
||||
import anndata
|
||||
import numpy as np
|
||||
import tiledb
|
||||
|
||||
from server.common.colors import convert_anndata_category_colors_to_cxg_category_colors
|
||||
from server.common.corpora import corpora_get_props_from_anndata
|
||||
from server.common.errors import ColorFormatException
|
||||
from server.common.utils.cxg_constants import CxgConstants
|
||||
from server.common.utils.cxg_generation_utils import (
|
||||
convert_dictionary_to_cxg_group,
|
||||
convert_dataframe_to_cxg_array,
|
||||
convert_ndarray_to_cxg_dense_array,
|
||||
convert_matrix_to_cxg_array,
|
||||
)
|
||||
from server.common.utils.matrix_utils import is_matrix_sparse, get_column_shift_encode_for_matrix
|
||||
|
||||
|
||||
class H5ADDataFile:
|
||||
""" Class encapsulating required information about an H5AD datafile that ultimately will be transformed into
|
||||
another format (currently just CXG is supported). """
|
||||
|
||||
def __init__(
|
||||
self,
|
||||
input_filename,
|
||||
backed=False,
|
||||
dataset_title=None,
|
||||
dataset_about=None,
|
||||
obs_index_column_name=None,
|
||||
vars_index_column_name=None,
|
||||
use_corpora_schema=True,
|
||||
):
|
||||
self.input_filename = input_filename
|
||||
self.backed = backed
|
||||
self.dataset_title = dataset_title
|
||||
self.dataset_about = dataset_about
|
||||
self.obs_index_column_name = obs_index_column_name
|
||||
self.vars_index_column_name = vars_index_column_name
|
||||
|
||||
self.use_corpora_schema = use_corpora_schema
|
||||
|
||||
self.validate_input_file_type()
|
||||
|
||||
self.extract_anndata_elements_from_file()
|
||||
self.extract_metadata_about_dataset()
|
||||
|
||||
self.validate_anndata()
|
||||
|
||||
def to_cxg(self, output_cxg_directory, sparse_threshold, convert_anndata_colors_to_cxg_colors=True):
|
||||
"""
|
||||
Writes the following attributes of the anndata to CXG: 1) the metadata as metadata attached to an empty
|
||||
DenseArray, 2) the obs DataFrame as a DenseArray, 3) the var DataFrame as a DenseArray, 4) all valid
|
||||
embeddings stored in obsm, each one as a DenseArray, 5) the main X matrix of the anndata as either a
|
||||
SparseArray or DenseArray based on the `sparse_threshold`, and optionally 6) the column shift of the main X
|
||||
matrix that might turn an otherwise Dense matrix into a Sparse matrix.
|
||||
"""
|
||||
|
||||
logging.info("Beginning writing to CXG.")
|
||||
ctx = tiledb.Ctx(
|
||||
{
|
||||
"sm.num_reader_threads": 32,
|
||||
"sm.num_writer_threads": 32,
|
||||
"sm.consolidation.buffer_size": 1 * 1024 * 1024 * 1024,
|
||||
}
|
||||
)
|
||||
|
||||
tiledb.group_create(output_cxg_directory, ctx=ctx)
|
||||
logging.info(f"\t...group created, with name {output_cxg_directory}")
|
||||
|
||||
convert_dictionary_to_cxg_group(
|
||||
output_cxg_directory, self.generate_cxg_metadata(convert_anndata_colors_to_cxg_colors)
|
||||
)
|
||||
logging.info("\t...dataset metadata saved")
|
||||
|
||||
convert_dataframe_to_cxg_array(output_cxg_directory, "obs", self.obs, self.obs_index_column_name, ctx)
|
||||
logging.info("\t...dataset obs dataframe saved")
|
||||
|
||||
convert_dataframe_to_cxg_array(output_cxg_directory, "var", self.var, self.var_index_column_name, ctx)
|
||||
logging.info("\t...dataset var dataframe saved")
|
||||
|
||||
self.write_anndata_embeddings_to_cxg(output_cxg_directory, ctx)
|
||||
logging.info("\t...dataset embeddings saved")
|
||||
|
||||
self.write_anndata_x_matrix_to_cxg(output_cxg_directory, ctx, sparse_threshold)
|
||||
logging.info("\t...dataset X matrix saved")
|
||||
|
||||
logging.info("Completed writing to CXG.")
|
||||
|
||||
def write_anndata_x_matrix_to_cxg(self, output_cxg_directory, ctx, sparse_threshold):
|
||||
matrix_container = f"{output_cxg_directory}/X"
|
||||
|
||||
x_matrix_data = self.anndata.X
|
||||
is_sparse = is_matrix_sparse(x_matrix_data, sparse_threshold)
|
||||
if not is_sparse:
|
||||
col_shift = get_column_shift_encode_for_matrix(x_matrix_data, sparse_threshold)
|
||||
is_sparse = col_shift is not None
|
||||
else:
|
||||
col_shift = None
|
||||
|
||||
if col_shift is not None:
|
||||
logging.info("Converting matrix X as sparse matrix with column shift encoding")
|
||||
x_col_shift_name = f"{output_cxg_directory}/X_col_shift"
|
||||
convert_ndarray_to_cxg_dense_array(x_col_shift_name, col_shift, ctx)
|
||||
|
||||
convert_matrix_to_cxg_array(matrix_container, x_matrix_data, is_sparse, ctx, col_shift)
|
||||
|
||||
tiledb.consolidate(matrix_container, ctx=ctx)
|
||||
if hasattr(tiledb, "vacuum"):
|
||||
tiledb.vacuum(matrix_container)
|
||||
|
||||
def write_anndata_embeddings_to_cxg(self, output_cxg_directory, ctx):
|
||||
def is_valid_embedding(adata, embedding_name, embedding_array):
|
||||
"""
|
||||
Returns true if this layout data is a valid array for front-end presentation with the following criteria:
|
||||
* ndarray, with shape (n_obs, >= 2), dtype float/int/uint
|
||||
* follows ScanPy embedding naming conventions
|
||||
* with all values finite or NaN (no +Inf or -Inf)
|
||||
"""
|
||||
|
||||
is_valid = isinstance(embedding_name, str) and embedding_name.startswith("X_") and len(embedding_name) > 2
|
||||
is_valid = is_valid and isinstance(embedding_array, np.ndarray) and embedding_array.dtype.kind in "fiu"
|
||||
is_valid = is_valid and embedding_array.shape[0] == adata.n_obs and embedding_array.shape[1] >= 2
|
||||
is_valid = is_valid and not np.any(np.isinf(embedding_array)) and not np.all(np.isnan(embedding_array))
|
||||
return is_valid
|
||||
|
||||
embedding_container = f"{output_cxg_directory}/emb"
|
||||
tiledb.group_create(embedding_container, ctx=ctx)
|
||||
|
||||
for embedding_name, embedding_values in self.anndata.obsm.items():
|
||||
if is_valid_embedding(self.anndata, embedding_name, embedding_values):
|
||||
embedding_name = f"{embedding_container}/{embedding_name[2:]}"
|
||||
convert_ndarray_to_cxg_dense_array(embedding_name, embedding_values, ctx)
|
||||
logging.info(f"\t\t...{embedding_name} embedding created")
|
||||
|
||||
def generate_cxg_metadata(self, convert_anndata_colors_to_cxg_colors):
|
||||
"""
|
||||
Return a dictionary containing metadata about CXG dataset. This include data about the version as well as
|
||||
Corpora schema properties if they exist, among other pieces of metadata.
|
||||
"""
|
||||
|
||||
cxg_group_metadata = {
|
||||
"cxg_version": CxgConstants.CXG_VERSION,
|
||||
"cxg_properties": json.dumps({"title": self.dataset_title, "about": self.dataset_about}),
|
||||
}
|
||||
if self.corpora_properties is not None:
|
||||
cxg_group_metadata["corpora"] = json.dumps(self.corpora_properties)
|
||||
|
||||
if convert_anndata_colors_to_cxg_colors:
|
||||
try:
|
||||
cxg_group_metadata["cxg_category_colors"] = json.dumps(
|
||||
convert_anndata_category_colors_to_cxg_category_colors(self.anndata)
|
||||
)
|
||||
except ColorFormatException:
|
||||
logging.warning(
|
||||
"Failed to extract colors from H5AD file! Fix the H5AD file or rerun with "
|
||||
"--disable-custom-colors. See help for more details."
|
||||
)
|
||||
|
||||
return cxg_group_metadata
|
||||
|
||||
def validate_input_file_type(self):
|
||||
"""
|
||||
Validate that the input file is of a type that we can handle. Currently the only valid file type is `.h5ad`.
|
||||
"""
|
||||
|
||||
if not self.input_filename.endswith(".h5ad"):
|
||||
raise Exception(f"Cannot process input file {self.input_filename}. File must be an H5AD.")
|
||||
|
||||
if self.dataset_title or self.dataset_about:
|
||||
logging.warning(
|
||||
"If you convert this dataset into CXG and you explicit specify values for the dataset title metadata "
|
||||
"or the dataset about metadata, it will override any metadata that is extracted as part of the "
|
||||
"Corpora schema fields."
|
||||
)
|
||||
|
||||
def validate_anndata(self):
|
||||
if not self.var.index.is_unique:
|
||||
raise ValueError("Variable index in AnnData object is not unique.")
|
||||
if not self.obs.index.is_unique:
|
||||
raise ValueError("Observation index in AnnData object is not unique.")
|
||||
|
||||
def extract_anndata_elements_from_file(self):
|
||||
logging.info(f"Reading in AnnData dataset: {path.basename(self.input_filename)}")
|
||||
self.anndata = anndata.read_h5ad(self.input_filename, backed="r" if self.backed else None)
|
||||
logging.info("Completed reading in AnnData dataset!")
|
||||
|
||||
self.obs = self.transform_dataframe_index_into_column(self.anndata.obs, "obs", self.obs_index_column_name)
|
||||
self.var = self.transform_dataframe_index_into_column(self.anndata.var, "var", self.vars_index_column_name)
|
||||
|
||||
def extract_metadata_about_dataset(self):
|
||||
"""
|
||||
Extract metadata information about the dataset that upon conversion will be saved as group metadata with the
|
||||
CXG that is generated. This metadata information includes Corpora schema properties, the dataset title and
|
||||
a link that details more information about the dataset.
|
||||
"""
|
||||
|
||||
self.corpora_properties = corpora_get_props_from_anndata(self.anndata) if self.use_corpora_schema else None
|
||||
if self.corpora_properties is None and self.use_corpora_schema:
|
||||
# If the return value is None, this means that we were not able to figure out what version of the Corpora
|
||||
# schema the object is using and therefore cannot extract any properties.
|
||||
raise ValueError("Unknown source file schema version is unsupported.")
|
||||
|
||||
# The title and about properties of the dataset are set by the following order: if they are explicitly defined
|
||||
# then use the explicit value. If the dataset is a Corpora-schema based schema, then extract the title and about
|
||||
# from the corpora_properties. Otherwise, use the input filename (only for title, about will be blank).
|
||||
if self.corpora_properties:
|
||||
corpora_project_links = self.corpora_properties.get("project_links", [])
|
||||
corpora_about_link = next(
|
||||
(link for link in corpora_project_links if (link.get("link_type", None) == "SUMMARY")), {}
|
||||
)
|
||||
else:
|
||||
corpora_about_link = {}
|
||||
|
||||
filename = path.splitext(path.basename(self.input_filename))[0]
|
||||
|
||||
self.dataset_title = self.dataset_title if self.dataset_title else corpora_about_link.get("link_name", filename)
|
||||
self.dataset_about = self.dataset_about if self.dataset_about else corpora_about_link.get("link_url")
|
||||
|
||||
def transform_dataframe_index_into_column(self, dataframe, dataframe_name, index_column_name):
|
||||
"""
|
||||
Convert the dataframe's index into another column in the dataframe. If an index_column_name is specified,
|
||||
use that column as the index instead.
|
||||
"""
|
||||
|
||||
if index_column_name is None:
|
||||
# Create a unique column name for the index.
|
||||
suffix = 0
|
||||
while f"name_{suffix}" in dataframe.columns:
|
||||
suffix += 1
|
||||
index_column_name = f"name_{suffix}"
|
||||
|
||||
# Turn the index into a normal column
|
||||
dataframe.rename_axis(index_column_name, inplace=True)
|
||||
dataframe.reset_index(inplace=True)
|
||||
|
||||
elif index_column_name in dataframe.columns:
|
||||
# User has specified alternative column for unique names, and it exists
|
||||
if not dataframe[index_column_name].is_unique:
|
||||
raise KeyError(
|
||||
f"Values in {dataframe_name}.{index_column_name} must be unique. Please prepare data to contain "
|
||||
f"unique values."
|
||||
)
|
||||
else:
|
||||
raise KeyError(f"Column {index_column_name} does not exist.")
|
||||
|
||||
setattr(self, f"{dataframe_name}_index_column_name", index_column_name)
|
||||
return dataframe
|
||||
@@ -423,11 +423,7 @@ class EndPointsAnndataAnnotations(unittest.TestCase, EndPointsAnnotations):
|
||||
cls.data, cls.tmp_dir, cls.annotations = data_with_tmp_annotations(
|
||||
MatrixDataType.H5AD, annotations_fixture=True
|
||||
)
|
||||
cls._setupClass(cls, [
|
||||
"--annotations-file",
|
||||
cls.annotations.output_file,
|
||||
cls.data.get_location(),
|
||||
])
|
||||
cls._setupClass(cls, ["--annotations-file", cls.annotations.output_file, cls.data.get_location(), ])
|
||||
|
||||
@classmethod
|
||||
def tearDownClass(cls):
|
||||
|
||||
@@ -1,9 +1,10 @@
|
||||
import unittest
|
||||
import anndata
|
||||
import json
|
||||
import tempfile
|
||||
import shutil
|
||||
import tempfile
|
||||
import unittest
|
||||
from http import HTTPStatus
|
||||
|
||||
import anndata
|
||||
import requests
|
||||
|
||||
from server.common.corpora import (
|
||||
@@ -104,7 +105,7 @@ class CorporaRESTAPITest(unittest.TestCase):
|
||||
"project_links": json.dumps([
|
||||
{"link_name": "test link", "link_type": "SUMMARY", "link_url": "https://a.u.r.l/"}
|
||||
]),
|
||||
"default_embedding": "X_tsne"
|
||||
"default_embedding": "X_tsne",
|
||||
}
|
||||
adata.uns.update(corpora_props)
|
||||
adata.write(path)
|
||||
|
||||
@@ -0,0 +1,149 @@
|
||||
import json
|
||||
import unittest
|
||||
from os import popen, path, mkdir
|
||||
from shutil import rmtree
|
||||
from uuid import uuid4
|
||||
|
||||
import numpy as np
|
||||
import tiledb
|
||||
from pandas import Series, DataFrame
|
||||
|
||||
from server.common.utils.cxg_generation_utils import (convert_dictionary_to_cxg_group, convert_dataframe_to_cxg_array,
|
||||
convert_ndarray_to_cxg_dense_array, convert_matrix_to_cxg_array)
|
||||
|
||||
PROJECT_ROOT = popen("git rev-parse --show-toplevel").read().strip()
|
||||
|
||||
|
||||
class TestCxgGenerationUtils(unittest.TestCase):
|
||||
def setUp(self):
|
||||
self.testing_cxg_temp_directory = f"{PROJECT_ROOT}/server/test/fixtures/{uuid4()}"
|
||||
mkdir(self.testing_cxg_temp_directory)
|
||||
|
||||
def tearDown(self):
|
||||
if path.isdir(self.testing_cxg_temp_directory):
|
||||
rmtree(self.testing_cxg_temp_directory)
|
||||
|
||||
def test__convert_dictionary_to_cxg_group__writes_successfully(self):
|
||||
random_dictionary = {"cookies": "chocolate_chip", "brownies": "chocolate", "cake": "double chocolate"}
|
||||
dictionary_name = "favorite_desserts"
|
||||
expected_array_directory = f"{self.testing_cxg_temp_directory}/{dictionary_name}"
|
||||
|
||||
convert_dictionary_to_cxg_group(self.testing_cxg_temp_directory, random_dictionary,
|
||||
group_metadata_name=dictionary_name)
|
||||
|
||||
array = tiledb.open(expected_array_directory)
|
||||
actual_stored_metadata = dict(array.meta.items())
|
||||
|
||||
self.assertTrue(path.isdir(expected_array_directory))
|
||||
self.assertTrue(isinstance(array, tiledb.DenseArray))
|
||||
self.assertEqual(random_dictionary, actual_stored_metadata)
|
||||
|
||||
def test__convert_dataframe_to_cxg_array__writes_successfully(self):
|
||||
random_int_category = Series(data=[3, 1, 2, 4], dtype=np.int64)
|
||||
random_bool_category = Series(data=[True, True, False, True], dtype=np.bool_)
|
||||
random_dataframe_name = f"random_dataframe_{uuid4()}"
|
||||
random_dataframe = DataFrame(data={"int_category": random_int_category, "bool_category": random_bool_category})
|
||||
|
||||
convert_dataframe_to_cxg_array(self.testing_cxg_temp_directory, random_dataframe_name, random_dataframe,
|
||||
"int_category", tiledb.Ctx())
|
||||
|
||||
expected_array_directory = f"{self.testing_cxg_temp_directory}/{random_dataframe_name}"
|
||||
expected_array_metadata = {
|
||||
"cxg_schema": json.dumps({"int_category": {"type": "int32"}, "bool_category": {"type": "boolean"},
|
||||
"index": "int_category"})}
|
||||
|
||||
actual_stored_dataframe_array = tiledb.open(expected_array_directory)
|
||||
actual_stored_dataframe_metadata = dict(actual_stored_dataframe_array.meta.items())
|
||||
|
||||
self.assertTrue(path.isdir(expected_array_directory))
|
||||
self.assertTrue(isinstance(actual_stored_dataframe_array, tiledb.DenseArray))
|
||||
self.assertDictEqual(expected_array_metadata, actual_stored_dataframe_metadata)
|
||||
self.assertTrue((actual_stored_dataframe_array[0:4]["int_category"] == random_int_category.to_numpy()).all())
|
||||
self.assertTrue((actual_stored_dataframe_array[0:4]["bool_category"] == random_bool_category.to_numpy()).all())
|
||||
|
||||
def test__convert_ndarray_to_cxg_dense_array__writes_successfully(self):
|
||||
ndarray = np.random.rand(3, 2)
|
||||
ndarray_name = f"{self.testing_cxg_temp_directory}/awesome_ndarray_{uuid4()}"
|
||||
|
||||
convert_ndarray_to_cxg_dense_array(ndarray_name, ndarray, tiledb.Ctx())
|
||||
|
||||
actual_stored_array = tiledb.open(ndarray_name)
|
||||
|
||||
self.assertTrue(path.isdir(ndarray_name))
|
||||
self.assertTrue(isinstance(actual_stored_array, tiledb.DenseArray))
|
||||
self.assertTrue((actual_stored_array[:, :] == ndarray).all())
|
||||
|
||||
def test__convert_matrix_to_cxg_array__dense_array_writes_successfully(self):
|
||||
matrix = np.float32(np.random.rand(3, 2))
|
||||
matrix_name = f"{self.testing_cxg_temp_directory}/awesome_matrix_{uuid4()}"
|
||||
|
||||
convert_matrix_to_cxg_array(matrix_name, matrix, False, tiledb.Ctx())
|
||||
|
||||
actual_stored_array = tiledb.open(matrix_name)
|
||||
|
||||
self.assertTrue(path.isdir(matrix_name))
|
||||
self.assertTrue(isinstance(actual_stored_array, tiledb.DenseArray))
|
||||
self.assertTrue((actual_stored_array[:, :] == matrix).all())
|
||||
|
||||
def test__convert_matrix_to_cxg_array__sparse_array_only_store_nonzeros_empty_array(self):
|
||||
matrix = np.zeros([3, 2])
|
||||
matrix_name = f"{self.testing_cxg_temp_directory}/awesome_zero_matrix_{uuid4()}"
|
||||
|
||||
convert_matrix_to_cxg_array(matrix_name, matrix, True, tiledb.Ctx())
|
||||
|
||||
actual_stored_array = tiledb.open(matrix_name)
|
||||
|
||||
self.assertTrue(path.isdir(matrix_name))
|
||||
self.assertTrue(isinstance(actual_stored_array, tiledb.SparseArray))
|
||||
self.assertTrue(actual_stored_array[:, :][''].size == 0)
|
||||
|
||||
def test__convert_matrix_to_cxg_array__sparse_array_only_store_nonzeros(self):
|
||||
matrix = np.zeros([3, 3])
|
||||
matrix[0, 0] = 1
|
||||
matrix[1, 1] = 1
|
||||
matrix[2, 2] = 2
|
||||
matrix_name = f"{self.testing_cxg_temp_directory}/awesome_sparse_matrix_{uuid4()}"
|
||||
|
||||
convert_matrix_to_cxg_array(matrix_name, matrix, True, tiledb.Ctx())
|
||||
|
||||
actual_stored_array = tiledb.open(matrix_name)
|
||||
|
||||
self.assertTrue(path.isdir(matrix_name))
|
||||
self.assertTrue(isinstance(actual_stored_array, tiledb.SparseArray))
|
||||
self.assertTrue(actual_stored_array[0, 0][''] == 1)
|
||||
self.assertTrue(actual_stored_array[1, 1][''] == 1)
|
||||
self.assertTrue(actual_stored_array[2, 2][''] == 2)
|
||||
self.assertTrue(actual_stored_array[:, :][''].size == 3)
|
||||
|
||||
def test__convert_matrix_to_cxg_array__sparse_array_with_column_encoding_empty_array(self):
|
||||
matrix_name = f"{self.testing_cxg_temp_directory}/awesome_column_shift_matrix_{uuid4()}"
|
||||
matrix = np.ones((3, 2))
|
||||
# The column shift will be equal to the matrix since subtracting the column shift from the matrix will create
|
||||
# a matrix of zeros which is sparse.
|
||||
column_shift = np.ones((3, 2))
|
||||
|
||||
convert_matrix_to_cxg_array(matrix_name, matrix, True, tiledb.Ctx(),
|
||||
column_shift_for_sparse_encoding=column_shift)
|
||||
|
||||
actual_stored_array = tiledb.open(matrix_name)
|
||||
|
||||
self.assertTrue(path.isdir(matrix_name))
|
||||
self.assertTrue(isinstance(actual_stored_array, tiledb.SparseArray))
|
||||
self.assertTrue(actual_stored_array[:, :][''].size == 0)
|
||||
|
||||
def test__convert_matrix_to_cxg_array__sparse_array_with_column_encoding_partial_array(self):
|
||||
matrix_name = f"{self.testing_cxg_temp_directory}/awesome_column_shift_matrix_{uuid4()}"
|
||||
matrix = np.ones((2, 2))
|
||||
# Only column shift the first column of ones.
|
||||
column_shift = np.array([[1, 0], [1, 0]])
|
||||
|
||||
convert_matrix_to_cxg_array(matrix_name, matrix, True, tiledb.Ctx(),
|
||||
column_shift_for_sparse_encoding=column_shift)
|
||||
|
||||
actual_stored_array = tiledb.open(matrix_name)
|
||||
|
||||
self.assertTrue(path.isdir(matrix_name))
|
||||
self.assertTrue(isinstance(actual_stored_array, tiledb.SparseArray))
|
||||
self.assertTrue(actual_stored_array[0, 1][''] == 1)
|
||||
self.assertTrue(actual_stored_array[1, 1][''] == 1)
|
||||
self.assertTrue(actual_stored_array[:, :][''].size == 2)
|
||||
@@ -2,10 +2,10 @@ import unittest
|
||||
from unittest.mock import patch
|
||||
|
||||
import numpy as np
|
||||
from pandas import Series
|
||||
from pandas import Series, DataFrame
|
||||
|
||||
from server.common.utils.type_conversion_utils import can_cast_to_float32, can_cast_to_int32, get_dtype_of_array, \
|
||||
get_schema_type_hint_of_array
|
||||
get_schema_type_hint_of_array, get_dtypes_and_schemas_of_dataframe
|
||||
|
||||
|
||||
class TestTypeConversionUtils(unittest.TestCase):
|
||||
@@ -119,3 +119,17 @@ class TestTypeConversionUtils(unittest.TestCase):
|
||||
i=test_type_index):
|
||||
array = Series(data=[], dtype=types[test_type_index])
|
||||
self.assertEqual(get_schema_type_hint_of_array(array), expected_schema_hints[test_type_index])
|
||||
|
||||
def test__get_dtypes_and_schemas_of_dataframe__dtype_and_schema_returns_as_expected(self):
|
||||
float_array = Series(data=[1, 2, 3], dtype=np.dtype(np.float64))
|
||||
category_array = Series(data=["a", "b", "b"], dtype="category")
|
||||
dataframe = DataFrame({"float_array": float_array, "category_array": category_array})
|
||||
|
||||
expected_data_types_dict = {"float_array": np.float32, "category_array": np.unicode}
|
||||
expected_schema_type_hints_dict = {"float_array": {"type": "float32"},
|
||||
"category_array": {"type": "categorical", "categories": ["a", "b"]}}
|
||||
|
||||
actual_dataframe_data_types, actual_dataframe_schema_type_hints = get_dtypes_and_schemas_of_dataframe(dataframe)
|
||||
|
||||
self.assertEqual(expected_data_types_dict, actual_dataframe_data_types)
|
||||
self.assertEqual(expected_schema_type_hints_dict, actual_dataframe_schema_type_hints)
|
||||
|
||||
@@ -1,14 +1,16 @@
|
||||
import os
|
||||
import tempfile
|
||||
import unittest
|
||||
from server.data_common.matrix_loader import MatrixDataLoader
|
||||
from server.test import PROJECT_ROOT, app_config, FIXTURES_ROOT
|
||||
|
||||
import numpy as np
|
||||
|
||||
import server.compute.diffexp_cxg as diffexp_cxg
|
||||
import server.compute.diffexp_generic as diffexp_generic
|
||||
from server.converters.cxgtool import write_cxg, create_cxg_group_metadata
|
||||
from server.test.performance.create_test_matrix import create_test_h5ad
|
||||
from server.converters.h5ad_data_file import H5ADDataFile
|
||||
from server.data_common.fbs.matrix import encode_matrix_fbs, decode_matrix_fbs
|
||||
import numpy as np
|
||||
import tempfile
|
||||
import os
|
||||
from server.data_common.matrix_loader import MatrixDataLoader
|
||||
from server.test import PROJECT_ROOT, app_config, FIXTURES_ROOT
|
||||
from server.test.performance.create_test_matrix import create_test_h5ad
|
||||
|
||||
|
||||
class DiffExpTest(unittest.TestCase):
|
||||
@@ -98,21 +100,22 @@ class DiffExpTest(unittest.TestCase):
|
||||
def sparse_diffexp(self, apply_col_shift):
|
||||
with tempfile.TemporaryDirectory() as dirname:
|
||||
# create a sparse matrix
|
||||
h5adfile = os.path.join(dirname, "sparse.h5ad")
|
||||
create_test_h5ad(h5adfile, 2000, 2000, 10, apply_col_shift)
|
||||
adaptor_anndata = self.load_dataset(h5adfile, extra_dataset_config=dict(embeddings__names=[]))
|
||||
adata = adaptor_anndata.data
|
||||
h5adfile_path = os.path.join(dirname, "sparse.h5ad")
|
||||
create_test_h5ad(h5adfile_path, 2000, 2000, 10, apply_col_shift)
|
||||
|
||||
h5ad_file_to_convert = H5ADDataFile(h5adfile_path, use_corpora_schema=False)
|
||||
|
||||
sparsename = os.path.join(dirname, "sparse.cxg")
|
||||
cxg_group_metadata = create_cxg_group_metadata(adata=adata, basefname="sparse.h5ad", title="sparse",)
|
||||
write_cxg(adata=adata, container=sparsename, cxg_group_metadata=cxg_group_metadata, sparse_threshold=11)
|
||||
h5ad_file_to_convert.to_cxg(sparsename, 11, True)
|
||||
|
||||
adaptor_anndata = self.load_dataset(h5adfile_path, extra_dataset_config=dict(embeddings__names=[]))
|
||||
|
||||
adaptor_sparse = self.load_dataset(sparsename)
|
||||
assert adaptor_sparse.open_array("X").schema.sparse
|
||||
assert adaptor_sparse.has_array("X_col_shift") == apply_col_shift
|
||||
|
||||
densename = os.path.join(dirname, "dense.cxg")
|
||||
cxg_group_metadata = create_cxg_group_metadata(adata=adata, basefname="dense.h5ad", title="dense",)
|
||||
write_cxg(adata=adata, container=densename, cxg_group_metadata=cxg_group_metadata, sparse_threshold=0)
|
||||
h5ad_file_to_convert.to_cxg(densename, True, 0)
|
||||
adaptor_dense = self.load_dataset(densename)
|
||||
assert not adaptor_dense.open_array("X").schema.sparse
|
||||
assert not adaptor_dense.has_array("X_col_shift")
|
||||
|
||||
@@ -1,41 +0,0 @@
|
||||
import shutil
|
||||
import unittest
|
||||
|
||||
import anndata
|
||||
|
||||
from server.common.data_locator import DataLocator
|
||||
from server.converters.cxgtool import write_cxg, create_cxg_group_metadata
|
||||
from server.data_cxg.cxg_adaptor import CxgAdaptor
|
||||
from server.test import PROJECT_ROOT, app_config, random_string
|
||||
from server.test.fixtures.fixtures import pbmc3k_colors
|
||||
|
||||
|
||||
class TestCxgAdaptor(unittest.TestCase):
|
||||
def setUp(self) -> None:
|
||||
self.fixtures = []
|
||||
|
||||
def tearDown(self) -> None:
|
||||
try:
|
||||
for data_locator in self.fixtures:
|
||||
print("REMOVING ", data_locator)
|
||||
shutil.rmtree(data_locator)
|
||||
except FileNotFoundError:
|
||||
pass
|
||||
|
||||
def test_cxg_category_colors(self):
|
||||
data = self.convert_pbmc3k(extract_colors=True)
|
||||
self.assertEqual(data.get_colors(), pbmc3k_colors)
|
||||
data = self.convert_pbmc3k(extract_colors=False)
|
||||
self.assertEqual(data.get_colors(), {})
|
||||
|
||||
def convert_pbmc3k(self, **kwargs):
|
||||
rand_str = random_string(8)
|
||||
data_locator = f"/tmp/test_{rand_str}.cxg"
|
||||
self.fixtures.append(data_locator)
|
||||
source_h5ad = anndata.read_h5ad(f"{PROJECT_ROOT}/example-dataset/pbmc3k.h5ad")
|
||||
cxg_group_metadata = create_cxg_group_metadata(
|
||||
adata=source_h5ad, basefname="pbmc3k.h5ad", title="pbmc3k", **kwargs
|
||||
)
|
||||
write_cxg(adata=source_h5ad, container=data_locator, cxg_group_metadata=cxg_group_metadata)
|
||||
config = app_config(data_locator)
|
||||
return CxgAdaptor(DataLocator(data_locator), config)
|
||||
@@ -0,0 +1,235 @@
|
||||
import json
|
||||
import unittest
|
||||
from glob import glob
|
||||
from os import popen, remove, path
|
||||
from shutil import rmtree
|
||||
from uuid import uuid4
|
||||
|
||||
import anndata
|
||||
import numpy as np
|
||||
from pandas import Series, DataFrame
|
||||
|
||||
from server.common.utils.corpora_constants import CorporaConstants
|
||||
from server.converters.h5ad_data_file import H5ADDataFile
|
||||
|
||||
PROJECT_ROOT = popen("git rev-parse --show-toplevel").read().strip()
|
||||
|
||||
|
||||
class TestH5ADDataFile(unittest.TestCase):
|
||||
|
||||
def setUp(self):
|
||||
self.sample_anndata = self._create_sample_anndata_dataset()
|
||||
self.sample_h5ad_filename = self._write_anndata_to_file(self.sample_anndata)
|
||||
|
||||
self.sample_output_directory = path.splitext(self.sample_h5ad_filename)[0] + ".cxg"
|
||||
|
||||
def tearDown(self):
|
||||
if self.sample_h5ad_filename:
|
||||
remove(self.sample_h5ad_filename)
|
||||
|
||||
if path.isdir(self.sample_output_directory):
|
||||
rmtree(self.sample_output_directory)
|
||||
|
||||
def test__create_h5ad_data_file__non_h5ad_raises_exception(self):
|
||||
non_h5ad_filename = "my_fancy_dataset.csv"
|
||||
|
||||
with self.assertRaises(Exception) as exception_context:
|
||||
H5ADDataFile(non_h5ad_filename)
|
||||
|
||||
self.assertIn("File must be an H5AD", str(exception_context.exception))
|
||||
|
||||
def test__create_h5ad_data_file__assert_warning_outputted_if_dataset_title_or_about_given(self):
|
||||
with self.assertLogs(level="WARN") as logger:
|
||||
H5ADDataFile(self.sample_h5ad_filename, dataset_title="My Awesome Dataset",
|
||||
dataset_about="http://www.awesomedataset.com", use_corpora_schema=False)
|
||||
|
||||
self.assertIn("will override any metadata that is extracted", logger.output[0])
|
||||
|
||||
def test__create_h5ad_data_file__reads_anndata_successfully(self):
|
||||
h5ad_file = H5ADDataFile(self.sample_h5ad_filename, use_corpora_schema=False)
|
||||
|
||||
self.assertTrue((h5ad_file.anndata.X == self.sample_anndata.X).all())
|
||||
self.assertEqual(h5ad_file.anndata.obs.sort_index(inplace=True),
|
||||
self.sample_anndata.obs.sort_index(inplace=True))
|
||||
self.assertEqual(h5ad_file.anndata.var.sort_index(inplace=True),
|
||||
self.sample_anndata.var.sort_index(inplace=True))
|
||||
|
||||
for key in h5ad_file.anndata.obsm.keys():
|
||||
self.assertIn(key, self.sample_anndata.obsm.keys())
|
||||
self.assertTrue((h5ad_file.anndata.obsm[key] == self.sample_anndata.obsm[key]).all())
|
||||
|
||||
for key in self.sample_anndata.obsm.keys():
|
||||
self.assertIn(key, h5ad_file.anndata.obsm.keys())
|
||||
self.assertTrue((h5ad_file.anndata.obsm[key] == self.sample_anndata.obsm[key]).all())
|
||||
|
||||
def test__create_h5ad_data_file__copies_index_of_obs_and_var_to_column(self):
|
||||
h5ad_file = H5ADDataFile(self.sample_h5ad_filename, use_corpora_schema=False)
|
||||
|
||||
# The automatic name chosen for the index should be "name_0"
|
||||
self.assertNotIn("name_0", self.sample_anndata.obs.columns)
|
||||
self.assertIn("name_0", h5ad_file.obs.columns)
|
||||
|
||||
self.assertNotIn("name_0", self.sample_anndata.var.columns)
|
||||
self.assertIn("name_0", h5ad_file.var.columns)
|
||||
|
||||
def test__create_h5ad_data_file__no_copy_if_obs_and_var_index_names_specified(self):
|
||||
h5ad_file = H5ADDataFile(self.sample_h5ad_filename, use_corpora_schema=False,
|
||||
obs_index_column_name="float_category", vars_index_column_name="int_category")
|
||||
|
||||
self.assertNotIn("name_0", h5ad_file.obs.columns)
|
||||
self.assertNotIn("name_0", h5ad_file.var.columns)
|
||||
|
||||
def test__create_h5ad_data_file__obs_and_var_index_names_specified_not_unique_raises_exception(self):
|
||||
|
||||
with self.assertRaises(Exception) as exception_context:
|
||||
H5ADDataFile(self.sample_h5ad_filename, use_corpora_schema=False,
|
||||
obs_index_column_name="float_category", vars_index_column_name="bool_category")
|
||||
|
||||
self.assertIn("Please prepare data to contain unique values", str(exception_context.exception))
|
||||
|
||||
def test__create_h5ad_data_file__obs_and_var_index_names_specified_doesnt_exist_raises_exception(self):
|
||||
with self.assertRaises(Exception) as exception_context:
|
||||
H5ADDataFile(self.sample_h5ad_filename, use_corpora_schema=False,
|
||||
obs_index_column_name="unknown_category", vars_index_column_name="i_dont_exist")
|
||||
|
||||
self.assertIn("does not exist", str(exception_context.exception))
|
||||
|
||||
def test__create_h5ad_data_file__extract_about_and_title_from_dataset(self):
|
||||
h5ad_file = H5ADDataFile(self.sample_h5ad_filename)
|
||||
|
||||
self.assertEqual(h5ad_file.dataset_title, "random_link_name")
|
||||
self.assertEqual(h5ad_file.dataset_about, "www.link.com")
|
||||
|
||||
def test__create_h5ad_data_file__inputted_dataset_title_and_about_overrides_extracted(self):
|
||||
h5ad_file = H5ADDataFile(self.sample_h5ad_filename, dataset_about="override_about",
|
||||
dataset_title="override_title")
|
||||
|
||||
self.assertEqual(h5ad_file.dataset_title, "override_title")
|
||||
self.assertEqual(h5ad_file.dataset_about, "override_about")
|
||||
|
||||
def test__to_cxg__simple_anndata_no_corpora_and_sparse(self):
|
||||
h5ad_file = H5ADDataFile(self.sample_h5ad_filename, use_corpora_schema=False)
|
||||
h5ad_file.to_cxg(self.sample_output_directory, 100)
|
||||
|
||||
self._validate_expected_generated_list_of_tiledb_files()
|
||||
|
||||
def test__to_cxg__simple_anndata_with_corpora_and_sparse(self):
|
||||
h5ad_file = H5ADDataFile(self.sample_h5ad_filename)
|
||||
h5ad_file.to_cxg(self.sample_output_directory, 100)
|
||||
|
||||
self._validate_expected_generated_list_of_tiledb_files()
|
||||
|
||||
def test__to_cxg__simple_anndata_no_corpora_and_dense(self):
|
||||
h5ad_file = H5ADDataFile(self.sample_h5ad_filename, use_corpora_schema=False)
|
||||
h5ad_file.to_cxg(self.sample_output_directory, 0)
|
||||
|
||||
self._validate_expected_generated_list_of_tiledb_files()
|
||||
|
||||
def test__to_cxg__simple_anndata_with_corpora_and_dense(self):
|
||||
h5ad_file = H5ADDataFile(self.sample_h5ad_filename)
|
||||
h5ad_file.to_cxg(self.sample_output_directory, 0)
|
||||
|
||||
self._validate_expected_generated_list_of_tiledb_files()
|
||||
|
||||
def test__to_cxg__with_sparse_column_encoding(self):
|
||||
anndata = self._create_sample_anndata_dataset()
|
||||
anndata.X = np.ones((3, 4))
|
||||
sparse_with_column_shift_filename = self._write_anndata_to_file(anndata)
|
||||
|
||||
h5ad_file = H5ADDataFile(sparse_with_column_shift_filename)
|
||||
h5ad_file.to_cxg(self.sample_output_directory, 50)
|
||||
|
||||
self._validate_expected_generated_list_of_tiledb_files(has_column_encoding=True)
|
||||
|
||||
# Clean up
|
||||
remove(sparse_with_column_shift_filename)
|
||||
|
||||
def _validate_expected_generated_list_of_tiledb_files(self, has_column_encoding=False):
|
||||
expected_directories, expected_obs_files, expected_var_files = \
|
||||
self._get_expected_generated_list_of_tiledb_files()
|
||||
|
||||
for directory in expected_directories:
|
||||
self.assertTrue(path.isdir(directory))
|
||||
|
||||
for obs_file in expected_obs_files:
|
||||
expected_location_of_obs_file = f"{self.sample_output_directory}/obs/*/{obs_file}"
|
||||
self.assertTrue(path.isfile(glob(expected_location_of_obs_file)[0]))
|
||||
|
||||
for var_file in expected_var_files:
|
||||
expected_location_of_var_file = f"{self.sample_output_directory}/var/*/{var_file}"
|
||||
self.assertTrue(path.isfile(glob(expected_location_of_var_file)[0]))
|
||||
|
||||
if has_column_encoding:
|
||||
self.assertTrue(path.isdir(f"{self.sample_output_directory}/X_col_shift"))
|
||||
|
||||
def _get_expected_generated_list_of_tiledb_files(self):
|
||||
|
||||
# Expected directories
|
||||
metadata_directory = f"{self.sample_output_directory}/cxg_group_metadata"
|
||||
main_x_directory = f"{self.sample_output_directory}/X"
|
||||
overall_embedding_directory = f"{self.sample_output_directory}/emb"
|
||||
specific_embedding_directory = f"{self.sample_output_directory}/emb/awesome_embedding"
|
||||
obs_directory = f"{self.sample_output_directory}/obs"
|
||||
var_directory = f"{self.sample_output_directory}/var"
|
||||
|
||||
# Obs files
|
||||
obs_files = []
|
||||
obs_files.append("name_0.tdb")
|
||||
obs_files.append("name_0_var.tdb")
|
||||
obs_files.append("string_category.tdb")
|
||||
obs_files.append("string_category_var.tdb")
|
||||
obs_files.append("float_category.tdb")
|
||||
|
||||
# Var files
|
||||
var_files = []
|
||||
var_files.append("name_0.tdb")
|
||||
var_files.append("name_0_var.tdb")
|
||||
var_files.append("bool_category.tdb")
|
||||
var_files.append("int_category.tdb")
|
||||
|
||||
return [metadata_directory, main_x_directory, overall_embedding_directory, specific_embedding_directory,
|
||||
obs_directory, var_directory], obs_files, var_files
|
||||
|
||||
def _write_anndata_to_file(self, anndata):
|
||||
temporary_filename = f"{PROJECT_ROOT}/server/test/fixtures/{uuid4()}.h5ad"
|
||||
anndata.write(temporary_filename)
|
||||
|
||||
return temporary_filename
|
||||
|
||||
def _create_sample_anndata_dataset(self):
|
||||
# Create X
|
||||
X = np.random.rand(3, 4)
|
||||
|
||||
# Create obs
|
||||
random_string_category = Series(data=["a", "b", "b"], dtype="category")
|
||||
random_float_category = Series(data=[3.2, 1.1, 2.2], dtype=np.float32)
|
||||
obs_dataframe = DataFrame(
|
||||
data={"string_category": random_string_category, "float_category": random_float_category})
|
||||
obs = obs_dataframe
|
||||
|
||||
# Create vars
|
||||
random_int_category = Series(data=[3, 1, 2, 4], dtype=np.int32)
|
||||
random_bool_category = Series(data=[True, True, False, True], dtype=np.bool_)
|
||||
var_dataframe = DataFrame(data={"int_category": random_int_category, "bool_category": random_bool_category})
|
||||
var = var_dataframe
|
||||
|
||||
# Create embeddings
|
||||
random_embedding = np.random.rand(3, 2)
|
||||
obsm = {"X_awesome_embedding": random_embedding}
|
||||
|
||||
# Create uns corpora metadata
|
||||
uns = {}
|
||||
for metadata_field in CorporaConstants.REQUIRED_SIMPLE_METADATA_FIELDS:
|
||||
uns[metadata_field] = "random"
|
||||
|
||||
for metadata_field in CorporaConstants.REQUIRED_JSON_ENCODED_METADATA_FIELD:
|
||||
uns[metadata_field] = json.dumps({"random_key": "random_value"})
|
||||
|
||||
# Need to carefully set the corpora schema versions in order for tests to pass.
|
||||
uns["version"] = {"corpora_schema_version": "1.0.0", "corpora_encoding_version": "0.1.0"}
|
||||
|
||||
# Set project links to be a dictionary
|
||||
uns["project_links"] = json.dumps(
|
||||
[{"link_name": "random_link_name", "link_url": "www.link.com", "link_type": "SUMMARY"}])
|
||||
|
||||
return anndata.AnnData(X=X, obs=obs, var=var, obsm=obsm, uns=uns)
|
||||
Reference in New Issue
Block a user