From 1b8fb974b63c505089b6fa2d49b1c82a3a091f7b Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Daniel=20F=C3=BCrth?= Date: Wed, 20 Sep 2023 16:31:05 +0200 Subject: [PATCH] updated cytosol training code --- .gitignore | 3 + augment.py | 10 +-- predict_cytosol.py | 208 +++++++++++++++++++++++++++++++++++++++++++++ predict_nuclei.py | 3 - train_cytosol.py | 161 +++++++++++++++++++++++++++++++++++ 5 files changed, 377 insertions(+), 8 deletions(-) create mode 100644 predict_cytosol.py create mode 100644 train_cytosol.py diff --git a/.gitignore b/.gitignore index 0ebd687..4129d8e 100644 --- a/.gitignore +++ b/.gitignore @@ -26,6 +26,9 @@ /training_data/dapi/images/*.tif /training_data/dapi/masks/*.tif +/training_data/cytosol/images/*.tif +/training_data/cytosol/masks/*.tif + # produced vignettes vignettes/*.html vignettes/*.pdf diff --git a/augment.py b/augment.py index 1931ae2..b437775 100644 --- a/augment.py +++ b/augment.py @@ -2,19 +2,19 @@ import cv2 # Set the paths for the input folders -mask_org_folder = "./training_data/masks_nuclei" -images_org_folder = "./training_data/images_nuclei" +mask_org_folder = "./training_data/masks_cytosol" +images_org_folder = "./training_data/images_cytosol" # Set the paths for the output folders -mask_folder = "./training_data/dapi/masks" -images_folder = "./training_data/dapi/images" +mask_folder = "./training_data/cytosol/masks" +images_folder = "./training_data/cytosol/images" # Create the output folders if they don't exist os.makedirs(mask_folder, exist_ok=True) os.makedirs(images_folder, exist_ok=True) # Set the ROI size and step -roi_size = (512, 512) +roi_size = (256, 256) roi_step = 128 # Get the list of file names in the mask_org folder diff --git a/predict_cytosol.py b/predict_cytosol.py new file mode 100644 index 0000000..feafbf5 --- /dev/null +++ b/predict_cytosol.py @@ -0,0 +1,208 @@ +#Load previously saved model +from keras.models import load_model +from tensorflow.keras.metrics import MeanIoU + +model = load_model("tutorial118_mitochondria_25epochs.hdf5", compile=False) + +import tensorflow as tf +import datetime + +from tensorflow.keras.utils import normalize +import os +import cv2 +from PIL import Image +import numpy as np +from matplotlib import pyplot as plt +from sklearn.preprocessing import MinMaxScaler +from keras.optimizers import Adam + +from glob import glob +from tqdm import tqdm +from tifffile import imread +from csbdeep.utils import Path, normalize + +image_directory = '/training_data/cytosol/images/' +mask_directory = '/training_data/cytosol/masks/' + +X = sorted(glob('training_data/cytosol/images/*.tif')) +Y = sorted(glob('training_data/cytosol/masks/*.tif')) +assert all(Path(x).name==Path(y).name for x,y in zip(X,Y)) + +image_names = sorted(glob('training_data/cytosol/images/*.tif')) + +SIZE = 256 +num_images = len(X) + +image_names_subset = image_names[0:num_images] + + +images = [cv2.imread(img, cv2.IMREAD_GRAYSCALE | cv2.IMREAD_ANYDEPTH) for img in image_names_subset] + + +image_dataset = np.array(images) +image_dataset = np.expand_dims(image_dataset, axis = 3) + +mask_names = sorted(glob("training_data/cytosol/masks/*.tif")) +mask_names_subset = mask_names[0:num_images] +masks = [cv2.imread(mask, 0) for mask in mask_names_subset] + +thresholded_masks = [] + +for mask in mask_names_subset: + # Load the image + mask_image = cv2.imread(mask, 0) + + # Threshold the image + _, thresholded_mask = cv2.threshold(mask_image, 1, 255, cv2.THRESH_BINARY) + + thresholded_masks.append(thresholded_mask) + +mask_dataset = np.array(thresholded_masks) +mask_dataset = np.expand_dims(mask_dataset, axis = 3) + +print("Image data shape is: ", image_dataset.shape) +print("Mask data shape is: ", mask_dataset.shape) +print("Max pixel value in image is: ", image_dataset.max()) +print("Labels in the mask are : ", np.unique(mask_dataset)) + +#Normalize images +image_dataset = image_dataset /image_dataset.max() #Can also normalize or scale using MinMax scaler +#Do not normalize masks, just rescale to 0 to 1. +mask_dataset = mask_dataset /255. #PIxel values will be 0 or 1 + +from sklearn.model_selection import train_test_split +X_train, X_test, y_train, y_test = train_test_split(image_dataset, mask_dataset, test_size = 0.20, random_state = 42) + + + +cv2_mat = cv2.cvtColor(X_test[0,:,:,:]) + + +# Assuming X_test is your NumPy array containing a 64-bit floating-point image +monochrome_image = X_test[0, :, :, 0] # Extract the first monochrome image + +# Convert the 64-bit floating-point image to 8-bit unsigned integer (0 to 255 range) +monochrome_image_uint8 = np.clip(monochrome_image * 255.0, 0, 255).astype(np.uint8) + +# Convert the 8-bit image to a BGR image +cv2_mat = cv2.cvtColor(monochrome_image_uint8, cv2.COLOR_GRAY2BGR) + +# Assuming you have cv2_mat as your image +cv2.imshow("Image", cv2_mat) +cv2.waitKey(0) # Add this to wait for a key press and then close the window +cv2.destroyAllWindows() # Add this to close all OpenCV windows when you're done + + +imread() + +#IOU +y_pred=model.predict(X_test) + + +float32_array = y_pred[0,:,:,:] + +float32_array = (float32_array > 0.5)*1.00 + + +# Assuming you have a NumPy array called 'float32_array' +# Scale the float32 array to the range [0, 255] +scaled_array = ((float32_array - float32_array.min()) / (float32_array.max() - float32_array.min()) * 255).astype(np.uint8) + +# Create a monochrome OpenCV Mat +monochrome_mat = cv2.cvtColor(scaled_array, cv2.COLOR_GRAY2BGR) + +# Convert the BGR image to grayscale +monochrome_mat = cv2.cvtColor(monochrome_mat, cv2.COLOR_BGR2GRAY) + + +# Assuming you have cv2_mat as your image +cv2.imshow("Image2", monochrome_mat) +cv2.waitKey(0) # Add this to wait for a key press and then close the window +cv2.destroyAllWindows() # Add this to close all OpenCV windows when you're done + + + +images = cv2.imread('./training_data/images_cytosol/CMYC_MAX_4.tif', cv2.IMREAD_GRAYSCALE | cv2.IMREAD_ANYDEPTH) + + +image_dataset = np.array(images) + +image_dataset = image_dataset /image_dataset.max() #Can also normalize or scale using MinMax scaler + + +# Assuming you have a NumPy array called 'image_dataset' with shape (2808, 2688) +original_height, original_width = image_dataset.shape +tile_size = 256 + +# Calculate the number of tiles in both dimensions +num_tiles_height = original_height // tile_size +num_tiles_width = original_width // tile_size + +# Initialize the tiled image array +tiledImage = np.zeros((num_tiles_height * num_tiles_width, tile_size, tile_size, 1), dtype=np.float32) + +# Iterate over the original image and extract tiles +tile_index = 0 +for i in range(num_tiles_height): + for j in range(num_tiles_width): + # Extract a 256x256 tile from the original image + tile = image_dataset[i * tile_size: (i + 1) * tile_size, j * tile_size: (j + 1) * tile_size] + tile = np.expand_dims(tile, axis=-1) # Add a channel dimension + tiledImage[tile_index] = tile + tile_index += 1 + + + +y_pred=model.predict(tiledImage) + + +# Assuming you have a NumPy array called 'tiledImage' with shape (x, 256, 256, 1) +num_tiles, tile_height, tile_width, num_channels = y_pred.shape +tile_size = 256 + +# Calculate the dimensions of the original image +original_height = num_tiles_height * tile_height +original_width = num_tiles_width * tile_width + +# Initialize the original image array +reconstructed_image = np.zeros((original_height, original_width, num_channels), dtype=np.float32) + +# Iterate over the tiles and reconstruct the original image +tile_index = 0 +for i in range(num_tiles_height): + for j in range(num_tiles_width): + # Get the tile from 'tiledImage' + tile = y_pred[tile_index] + + # Place the tile in the reconstructed image + reconstructed_image[i * tile_size: (i + 1) * tile_size, j * tile_size: (j + 1) * tile_size] = tile + tile_index += 1 + + +float32_array = reconstructed_image + +# Assuming you have a NumPy array called 'float32_array' +# Scale the float32 array to the range [0, 255] +scaled_array = ((float32_array - float32_array.min()) / (float32_array.max() - float32_array.min()) * 255).astype(np.uint8) + +# Create a monochrome OpenCV Mat +monochrome_mat = cv2.cvtColor(scaled_array, cv2.COLOR_GRAY2BGR) + +# Convert the BGR image to grayscale +monochrome_mat = cv2.cvtColor(monochrome_mat, cv2.COLOR_BGR2GRAY) + + +# Assuming you have cv2_mat as your image +cv2.imshow("ImageNew", monochrome_mat) +cv2.waitKey(0) # Add this to wait for a key press and then close the window +cv2.destroyAllWindows() # Add this to close all OpenCV windows when you're done + + +# Assuming you have an OpenCV Mat called 'image_mat' that you want to save as a tiled TIFF +output_file_path = "output.tiff" # Specify the path where you want to save the tiled TIFF + + +success = cv2.imwrite(output_file_path, monochrome_mat) + + + diff --git a/predict_nuclei.py b/predict_nuclei.py index a81aa31..32e7526 100644 --- a/predict_nuclei.py +++ b/predict_nuclei.py @@ -13,9 +13,6 @@ from stardist import random_label_cmap, _draw_polygons, export_imagej_rois from stardist.models import StarDist2D -np.random.seed(6) -lbl_cmap = random_label_cmap() - X = sorted(glob('./training_data/images_nuclei/*.tif')) X = list(map(imread,X)) diff --git a/train_cytosol.py b/train_cytosol.py new file mode 100644 index 0000000..d51c454 --- /dev/null +++ b/train_cytosol.py @@ -0,0 +1,161 @@ +import tensorflow as tf +import datetime + +from tensorflow.keras.utils import normalize +import os +import cv2 +from PIL import Image +import numpy as np +from matplotlib import pyplot as plt +from sklearn.preprocessing import MinMaxScaler +from keras.optimizers import Adam + +from glob import glob +from tqdm import tqdm +from tifffile import imread +from csbdeep.utils import Path, normalize + +image_directory = '/training_data/cytosol/images/' +mask_directory = '/training_data/cytosol/masks/' + +X = sorted(glob('training_data/cytosol/images/*.tif')) +Y = sorted(glob('training_data/cytosol/masks/*.tif')) +assert all(Path(x).name==Path(y).name for x,y in zip(X,Y)) + +image_names = sorted(glob('training_data/cytosol/images/*.tif')) + +SIZE = 256 +num_images = len(X) + +image_names_subset = image_names[0:num_images] + + +images = [cv2.imread(img, cv2.IMREAD_GRAYSCALE | cv2.IMREAD_ANYDEPTH) for img in image_names_subset] + + +image_dataset = np.array(images) +image_dataset = np.expand_dims(image_dataset, axis = 3) + +mask_names = sorted(glob("training_data/cytosol/masks/*.tif")) +mask_names_subset = mask_names[0:num_images] +masks = [cv2.imread(mask, 0) for mask in mask_names_subset] + +thresholded_masks = [] + +for mask in mask_names_subset: + # Load the image + mask_image = cv2.imread(mask, 0) + + # Define the kernel (structuring element) for erosion + kernel_size = 3 # Adjust the size as needed + kernel = np.ones((kernel_size, kernel_size), dtype=np.uint8) + + # Perform erosion + eroded_image = cv2.erode(mask_image, kernel, iterations=3) + + # Threshold the image + _, thresholded_mask = cv2.threshold(eroded_image, 1, 255, cv2.THRESH_BINARY) + + thresholded_masks.append(thresholded_mask) + +mask_dataset = np.array(thresholded_masks) +mask_dataset = np.expand_dims(mask_dataset, axis = 3) + +print("Image data shape is: ", image_dataset.shape) +print("Mask data shape is: ", mask_dataset.shape) +print("Max pixel value in image is: ", image_dataset.max()) +print("Labels in the mask are : ", np.unique(mask_dataset)) + +#Normalize images +image_dataset = image_dataset /image_dataset.max() #Can also normalize or scale using MinMax scaler +#Do not normalize masks, just rescale to 0 to 1. +mask_dataset = mask_dataset /255. #PIxel values will be 0 or 1 + +from sklearn.model_selection import train_test_split +X_train, X_test, y_train, y_test = train_test_split(image_dataset, mask_dataset, test_size = 0.20, random_state = 42) + +# Building Unet by dividing encoder and decoder into blocks + +from keras.models import Model +from keras.layers import Input, Conv2D, MaxPooling2D, UpSampling2D, concatenate, Conv2DTranspose, BatchNormalization, Dropout, Lambda +from keras.optimizers import Adam +from keras.layers import Activation, MaxPool2D, Concatenate + + +def conv_block(input, num_filters): + x = Conv2D(num_filters, 3, padding="same")(input) + x = BatchNormalization()(x) #Not in the original network. + x = Activation("relu")(x) + + x = Conv2D(num_filters, 3, padding="same")(x) + x = BatchNormalization()(x) #Not in the original network + x = Activation("relu")(x) + + return x + +#Encoder block: Conv block followed by maxpooling + + +def encoder_block(input, num_filters): + x = conv_block(input, num_filters) + p = MaxPool2D((2, 2))(x) + return x, p + +#Decoder block +#skip features gets input from encoder for concatenation + +def decoder_block(input, skip_features, num_filters): + x = Conv2DTranspose(num_filters, (2, 2), strides=2, padding="same")(input) + x = Concatenate()([x, skip_features]) + x = conv_block(x, num_filters) + return x + +#Build Unet using the blocks +def build_unet(input_shape, n_classes): + inputs = Input(input_shape) + + s1, p1 = encoder_block(inputs, 64) + s2, p2 = encoder_block(p1, 128) + s3, p3 = encoder_block(p2, 256) + s4, p4 = encoder_block(p3, 512) + + b1 = conv_block(p4, 1024) #Bridge + + d1 = decoder_block(b1, s4, 512) + d2 = decoder_block(d1, s3, 256) + d3 = decoder_block(d2, s2, 128) + d4 = decoder_block(d3, s1, 64) + + if n_classes == 1: #Binary + activation = 'sigmoid' + else: + activation = 'softmax' + + outputs = Conv2D(n_classes, 1, padding="same", activation=activation)(d4) #Change the activation based on n_classes + print(activation) + + model = Model(inputs, outputs, name="U-Net") + return model + +IMG_HEIGHT = image_dataset.shape[1] +IMG_WIDTH = image_dataset.shape[2] +IMG_CHANNELS = image_dataset.shape[3] + +input_shape = (IMG_HEIGHT, IMG_WIDTH, IMG_CHANNELS) + +model = build_unet(input_shape, n_classes=1) +model.compile(optimizer=Adam(learning_rate = 1e-3), loss='binary_crossentropy', metrics=['accuracy']) +model.summary() + +log_dir = "logs/fit/" + datetime.datetime.now().strftime("%Y%m%d-%H%M%S") +tensorboard_callback = tf.keras.callbacks.TensorBoard(log_dir=log_dir, histogram_freq=1) + +history = model.fit(X_train, y_train, + batch_size = 8, + verbose=1, + epochs=100, + validation_data=(X_test, y_test), + callbacks=[tensorboard_callback], + shuffle=False) + +model.save('tutorial118_mitochondria_100epochs.hdf5')