update
Browse files
scDifformer_data_process/loom_batch_100.py
CHANGED
|
@@ -18,32 +18,32 @@ def set_log() -> None:
|
|
| 18 |
set_log()
|
| 19 |
|
| 20 |
|
| 21 |
-
#
|
| 22 |
input_folder = '/mnt/nfs/data/geneformer/output/v2/looms'
|
| 23 |
output_folder = '/mnt/nfs/data/geneformer/output/v2/looms_sub'
|
| 24 |
per_group_num = 100
|
| 25 |
|
| 26 |
|
| 27 |
if __name__ == '__main__':
|
| 28 |
-
#
|
| 29 |
os.makedirs(output_folder, exist_ok=True)
|
| 30 |
|
| 31 |
-
#
|
| 32 |
files = os.listdir(input_folder)
|
| 33 |
group_count = 0
|
| 34 |
|
| 35 |
for i, file in enumerate(files, start=1):
|
| 36 |
-
#
|
| 37 |
group_folder = os.path.join(output_folder, f'group_{group_count + 1}')
|
| 38 |
|
| 39 |
-
#
|
| 40 |
os.makedirs(group_folder, exist_ok=True)
|
| 41 |
|
| 42 |
-
#
|
| 43 |
input_file = os.path.join(input_folder, file)
|
| 44 |
output_file = os.path.join(group_folder, file)
|
| 45 |
|
| 46 |
-
#
|
| 47 |
shutil.copy(input_file, output_file)
|
| 48 |
|
| 49 |
if i % per_group_num == 0:
|
|
|
|
| 18 |
set_log()
|
| 19 |
|
| 20 |
|
| 21 |
+
# Input and output folder paths
|
| 22 |
input_folder = '/mnt/nfs/data/geneformer/output/v2/looms'
|
| 23 |
output_folder = '/mnt/nfs/data/geneformer/output/v2/looms_sub'
|
| 24 |
per_group_num = 100
|
| 25 |
|
| 26 |
|
| 27 |
if __name__ == '__main__':
|
| 28 |
+
# Create output folder if it does not exist
|
| 29 |
os.makedirs(output_folder, exist_ok=True)
|
| 30 |
|
| 31 |
+
# Iterate through files in the input folder
|
| 32 |
files = os.listdir(input_folder)
|
| 33 |
group_count = 0
|
| 34 |
|
| 35 |
for i, file in enumerate(files, start=1):
|
| 36 |
+
# Build output subfolder path
|
| 37 |
group_folder = os.path.join(output_folder, f'group_{group_count + 1}')
|
| 38 |
|
| 39 |
+
# Create output subfolder if it does not exist
|
| 40 |
os.makedirs(group_folder, exist_ok=True)
|
| 41 |
|
| 42 |
+
# Build input and output file paths
|
| 43 |
input_file = os.path.join(input_folder, file)
|
| 44 |
output_file = os.path.join(group_folder, file)
|
| 45 |
|
| 46 |
+
# Copy file into the output subfolder
|
| 47 |
shutil.copy(input_file, output_file)
|
| 48 |
|
| 49 |
if i % per_group_num == 0:
|
scDifformer_data_process/preflight/convert_data.py
CHANGED
|
@@ -25,7 +25,7 @@ gene_info_table_path = '/mnt/nfs/data/geneformer/input/material/gene_info_table.
|
|
| 25 |
|
| 26 |
gene_info = pd.read_csv(gene_info_table_path)
|
| 27 |
|
| 28 |
-
#
|
| 29 |
gene_name_id_combine_dict = gene_info.set_index('gene_name')['ensembl_id'].to_dict()
|
| 30 |
|
| 31 |
input_dir = '/mnt/nfs/hs/disco_1115_all'
|
|
|
|
| 25 |
|
| 26 |
gene_info = pd.read_csv(gene_info_table_path)
|
| 27 |
|
| 28 |
+
# Convert to dictionary
|
| 29 |
gene_name_id_combine_dict = gene_info.set_index('gene_name')['ensembl_id'].to_dict()
|
| 30 |
|
| 31 |
input_dir = '/mnt/nfs/hs/disco_1115_all'
|
scDifformer_data_process/preprocess/preprocess.py
CHANGED
|
@@ -32,7 +32,7 @@ home_dir_list = ['/scDifformer/data/240412_test/raw_h5ad']
|
|
| 32 |
|
| 33 |
gene_info = pd.read_csv(gene_info_table_path)
|
| 34 |
|
| 35 |
-
#
|
| 36 |
gene_name_id_combine_dict = gene_info.set_index('gene_name')['ensembl_id'].to_dict()
|
| 37 |
gene_name_type_dict = gene_info.set_index('gene_name')['gene_type'].to_dict()
|
| 38 |
gene_id_type_dict = gene_info.set_index('ensembl_id')['gene_type'].to_dict()
|
|
@@ -47,34 +47,34 @@ def preprocess(file_path: str) -> None:
|
|
| 47 |
adata = sc.read_h5ad(file_path)
|
| 48 |
adata.raw = None
|
| 49 |
adata.X = np.round(adata.X, decimals=1)
|
| 50 |
-
# 1.
|
| 51 |
adata.var['mt'] = adata.var_names.str.startswith('MT-') # annotate the group of mitochondrial genes as 'mt'
|
| 52 |
sc.pp.calculate_qc_metrics(adata, qc_vars=['mt'], percent_top=None, log1p=False, inplace=True)
|
| 53 |
|
| 54 |
-
#
|
| 55 |
mean_ratio_mt = np.mean(adata.obs["pct_counts_mt"])
|
| 56 |
std_ratio_mt = np.std(adata.obs["pct_counts_mt"])
|
| 57 |
|
| 58 |
-
#
|
| 59 |
min_ratio_mt = mean_ratio_mt - 3 * std_ratio_mt
|
| 60 |
max_ratio_mt = mean_ratio_mt + 3 * std_ratio_mt
|
| 61 |
within_range_mt = (adata.obs["pct_counts_mt"] >= min_ratio_mt) & (adata.obs["pct_counts_mt"] <= max_ratio_mt)
|
| 62 |
-
#
|
| 63 |
adata = adata[within_range_mt, :]
|
| 64 |
|
| 65 |
-
#
|
| 66 |
mean_ratio = np.mean(adata.obs["total_counts"])
|
| 67 |
std_ratio = np.std(adata.obs["total_counts"])
|
| 68 |
|
| 69 |
-
#
|
| 70 |
min_ratio = mean_ratio - 3 * std_ratio
|
| 71 |
max_ratio = mean_ratio + 3 * std_ratio
|
| 72 |
within_range = (adata.obs["total_counts"] >= min_ratio) & (adata.obs["total_counts"] <= max_ratio)
|
| 73 |
|
| 74 |
-
#
|
| 75 |
adata = adata[within_range]
|
| 76 |
|
| 77 |
-
# 2.
|
| 78 |
scrub = scrublet.Scrublet(adata.X)
|
| 79 |
# Run doublet detection
|
| 80 |
doublet_scores, predicted_doublets = scrub.scrub_doublets()
|
|
@@ -141,6 +141,6 @@ if __name__ == '__main__':
|
|
| 141 |
home_path = home_dir
|
| 142 |
files = [os.path.join(home_path, file_name) for file_name in os.listdir(home_path) if
|
| 143 |
str(file_name).endswith(".h5ad") and file_name not in no_need_handle_set]
|
| 144 |
-
logger.info(f'
|
| 145 |
with mp.Pool(processes=32) as pool:
|
| 146 |
pool.map(preprocess, files)
|
|
|
|
| 32 |
|
| 33 |
gene_info = pd.read_csv(gene_info_table_path)
|
| 34 |
|
| 35 |
+
# Convert to dictionaries
|
| 36 |
gene_name_id_combine_dict = gene_info.set_index('gene_name')['ensembl_id'].to_dict()
|
| 37 |
gene_name_type_dict = gene_info.set_index('gene_name')['gene_type'].to_dict()
|
| 38 |
gene_id_type_dict = gene_info.set_index('ensembl_id')['gene_type'].to_dict()
|
|
|
|
| 47 |
adata = sc.read_h5ad(file_path)
|
| 48 |
adata.raw = None
|
| 49 |
adata.X = np.round(adata.X, decimals=1)
|
| 50 |
+
# 1. Filter mitochondrial ratio and total gene counts within 3 standard deviations of the mean
|
| 51 |
adata.var['mt'] = adata.var_names.str.startswith('MT-') # annotate the group of mitochondrial genes as 'mt'
|
| 52 |
sc.pp.calculate_qc_metrics(adata, qc_vars=['mt'], percent_top=None, log1p=False, inplace=True)
|
| 53 |
|
| 54 |
+
# Compute mean and standard deviation of mitochondrial ratio
|
| 55 |
mean_ratio_mt = np.mean(adata.obs["pct_counts_mt"])
|
| 56 |
std_ratio_mt = np.std(adata.obs["pct_counts_mt"])
|
| 57 |
|
| 58 |
+
# Keep cells within three standard deviations of the mean
|
| 59 |
min_ratio_mt = mean_ratio_mt - 3 * std_ratio_mt
|
| 60 |
max_ratio_mt = mean_ratio_mt + 3 * std_ratio_mt
|
| 61 |
within_range_mt = (adata.obs["pct_counts_mt"] >= min_ratio_mt) & (adata.obs["pct_counts_mt"] <= max_ratio_mt)
|
| 62 |
+
# Remove cells outside the valid range
|
| 63 |
adata = adata[within_range_mt, :]
|
| 64 |
|
| 65 |
+
# Compute mean and standard deviation of total gene counts
|
| 66 |
mean_ratio = np.mean(adata.obs["total_counts"])
|
| 67 |
std_ratio = np.std(adata.obs["total_counts"])
|
| 68 |
|
| 69 |
+
# Keep cells within three standard deviations of the mean
|
| 70 |
min_ratio = mean_ratio - 3 * std_ratio
|
| 71 |
max_ratio = mean_ratio + 3 * std_ratio
|
| 72 |
within_range = (adata.obs["total_counts"] >= min_ratio) & (adata.obs["total_counts"] <= max_ratio)
|
| 73 |
|
| 74 |
+
# Remove cells outside the valid range
|
| 75 |
adata = adata[within_range]
|
| 76 |
|
| 77 |
+
# 2. Create a DoubletFinder-like object and run doublet detection
|
| 78 |
scrub = scrublet.Scrublet(adata.X)
|
| 79 |
# Run doublet detection
|
| 80 |
doublet_scores, predicted_doublets = scrub.scrub_doublets()
|
|
|
|
| 141 |
home_path = home_dir
|
| 142 |
files = [os.path.join(home_path, file_name) for file_name in os.listdir(home_path) if
|
| 143 |
str(file_name).endswith(".h5ad") and file_name not in no_need_handle_set]
|
| 144 |
+
logger.info(f'Number of files processed in this run: {len(files)}')
|
| 145 |
with mp.Pool(processes=32) as pool:
|
| 146 |
pool.map(preprocess, files)
|
scDifformer_data_process/print_data_info.py
CHANGED
|
@@ -12,22 +12,22 @@
|
|
| 12 |
#
|
| 13 |
# file_path = '/scDifformer/data/240412_test/looms/zheng68k_fold0_test.loom'
|
| 14 |
#
|
| 15 |
-
# #
|
| 16 |
# with loompy.connect(file_path) as ds:
|
| 17 |
#
|
| 18 |
-
# #
|
| 19 |
# print("\nColumn:")
|
| 20 |
# for key, val in ds.ca.items():
|
| 21 |
# print(f"{key}: {val}")
|
| 22 |
#
|
| 23 |
-
# #
|
| 24 |
# print("\nRow:")
|
| 25 |
# for key, val in ds.ra.items():
|
| 26 |
# print(f"{key}: {val}")
|
| 27 |
#
|
| 28 |
-
# #
|
| 29 |
# print("\nData matrix:")
|
| 30 |
-
# matrix_slice = ds[:, :5] #
|
| 31 |
# print(matrix_slice)
|
| 32 |
|
| 33 |
# import pickle
|
|
|
|
| 12 |
#
|
| 13 |
# file_path = '/scDifformer/data/240412_test/looms/zheng68k_fold0_test.loom'
|
| 14 |
#
|
| 15 |
+
# # Use loompy to connect to the loom file
|
| 16 |
# with loompy.connect(file_path) as ds:
|
| 17 |
#
|
| 18 |
+
# # Print all column attributes (cell metadata)
|
| 19 |
# print("\nColumn:")
|
| 20 |
# for key, val in ds.ca.items():
|
| 21 |
# print(f"{key}: {val}")
|
| 22 |
#
|
| 23 |
+
# # Print all row attributes (gene metadata)
|
| 24 |
# print("\nRow:")
|
| 25 |
# for key, val in ds.ra.items():
|
| 26 |
# print(f"{key}: {val}")
|
| 27 |
#
|
| 28 |
+
# # Print part of the matrix (e.g., expression of first 5 genes and first 5 cells)
|
| 29 |
# print("\nData matrix:")
|
| 30 |
+
# matrix_slice = ds[:, :5] # Expression of all genes in the first 5 cells
|
| 31 |
# print(matrix_slice)
|
| 32 |
|
| 33 |
# import pickle
|
scDifformer_data_process/tokenizer.py
CHANGED
|
@@ -127,7 +127,7 @@ class TranscriptomeTokenizer:
|
|
| 127 |
i = 1
|
| 128 |
# loops through directories to tokenize .loom files
|
| 129 |
for loom_file_path in loom_data_directory.glob("*.loom"):
|
| 130 |
-
logger.info(f"
|
| 131 |
i += 1
|
| 132 |
print(f"Tokenizing {loom_file_path}")
|
| 133 |
file_tokenized_cells, file_cell_metadata = self.tokenize_file(
|
|
|
|
| 127 |
i = 1
|
| 128 |
# loops through directories to tokenize .loom files
|
| 129 |
for loom_file_path in loom_data_directory.glob("*.loom"):
|
| 130 |
+
logger.info(f"Number of files processed: {i}")
|
| 131 |
i += 1
|
| 132 |
print(f"Tokenizing {loom_file_path}")
|
| 133 |
file_tokenized_cells, file_cell_metadata = self.tokenize_file(
|