<jupyter_start><jupyter_text>Predicting surface location using fluorescent beads Fluorescent beads on a clean coverslip were imaged with confocal laser scanning microscope. Their 3D position have been determined using 3D Gaussian fit. Given their radius, this is enough information to uncover coverslip surface position beneath them. However, bead radius has a wide variation. In this file, we learn how to correct for bead size using various parameters that 3D Gaussian fit spits out and an independent measurment of actual coverslip surface position.  Prepare data<jupyter_code># Use seaborn for pairplot
!pip install -q seaborn
# Use some functions from tensorflow_docs
!pip install -q git+https://github.com/tensorflow/docs
import pathlib
import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
import seaborn as sns
import tensorflow as tf
from tensorflow import keras
from tensorflow.keras import layers
import tensorflow_docs as tfdocs
import tensorflow_docs.plots
import tensorflow_docs.modeling<jupyter_output>Building wheel for tensorflow-docs (setup.py) ... [?25l[?25hdone<jupyter_text>Download and save data<jupyter_code>!git clone https://github.com/masterika/Fluorescence_BSVC.git
import os
# Root directory of the project
ROOT_DIR = os.getcwd()
os.chdir('Fluorescence_BSVC/Datasets')
!unzip coverslip.zip # Results in creating a folder named "bacteria"
bacteria_dataset_path = os.path.join(os.getcwd(), "bacteria")
training_set = pd.read_csv("training.csv")
test_set = pd.read_csv("test.csv")
os.chdir(ROOT_DIR) # go back to root directory<jupyter_output>Cloning into 'Fluorescence_BSVC'...
remote: Enumerating objects: 55, done.[K
remote: Counting objects: 100% (55/55), done.[K
remote: Compressing objects: 100% (50/50), done.[K
remote: Total 55 (delta 18), reused 0 (delta 0), pack-reused 0[K
Unpacking objects: 100% (55/55), done.
Archive:  coverslip.zip
  inflating: test.csv                
  inflating: training.csv<jupyter_text>Inspect data<jupyter_code>training_set.head()<jupyter_output><empty_output><jupyter_text>Column contents are as follows:- x, y, and z: three dimensional position of fluorescent bead's center, 3D Gaussian fit result- stdx, stdy and stdz: gaussian spread (standard deviation) of fluorescent bead intensity in x, y and z directions respectively, fit result- prop: one of the fit coefficients, overal intensity indicator- const: Offset, fit coefficient- volume: number of voxels comprising bead- area: number of pixels surrounding the bead- dif: the label. Relative z position of bead's center to coverslip surface obtained by substracting coverslip position in bead's vicinity from bead's z position. This uses an independent measurement of coverslip surface positionWith 3D Gaussian fit: $$ const+prop\times e^\frac{(x-x_{av})^2}{2stdx^2}e^\frac{(y-y_{av})^2}{2stdy^2}e^\frac{(z-z_{av})^2}{2stdz^2}$$<jupyter_code>sns.pairplot(training_set[["dif","stdz","stdx","volume","prop"]], diag_kind="kde")<jupyter_output><empty_output><jupyter_text>Separate labels from data<jupyter_code>LABEL = "dif"	
training_labels = training_set.pop(LABEL)
test_labels = test_set.pop(LABEL)<jupyter_output><empty_output><jupyter_text>Normalize data<jupyter_code>training_set_n=(training_set-training_set.mean())/training_set.std()
test_set_n=(test_set-training_set.mean())/training_set.std()
# labels are already pretty close to 0-1 range, so we don't normilize them.<jupyter_output><empty_output><jupyter_text>Models<jupyter_code>FEATURES = ["stdz","stdx","volume","prop","stdy","area","const"] # set of meaniningful features for dif prediction. stdz feature does most of the job, but rest help.<jupyter_output><empty_output><jupyter_text>Based on pair-plots above, dif seems to be approximately linear in features. Therefore, first let's see how the most intuitive, linear model will do<jupyter_code>def build_linear_model():
  model = keras.Sequential([        
    layers.Dense(1)
  ])

  optimizer = tf.keras.optimizers.Adam(0.001)

  model.compile(loss='mse',
                optimizer=optimizer,
                metrics=['mae'])
  return model<jupyter_output><empty_output><jupyter_text>Looking closer, plots seem consistently a bit curved. Let's see if neural networks can detect these subtle dependencies and improve precision<jupyter_code>def build_neural_network():
  model = keras.Sequential([
    layers.Dense(6, activation='relu'),    
    layers.Dense(1)
  ])

  optimizer = tf.keras.optimizers.Adam(0.001)

  model.compile(loss='mse',
                optimizer=optimizer,
                metrics=['mae'])
  return model<jupyter_output><empty_output><jupyter_text>Training with linear model<jupyter_code>model = build_linear_model() 
EPOCHS = 2000
early_stop = keras.callbacks.EarlyStopping(monitor='val_loss', patience=10) # stop training whenever validation set stops showing improvement
early_history = model.fit(training_set_n[FEATURES], training_labels, 
                    epochs=EPOCHS, validation_split = 0.3, verbose=0, 
                    callbacks=[early_stop, tfdocs.modeling.EpochDots()])<jupyter_output>WARNING:tensorflow:Layer dense_12 is casting an input tensor from dtype float64 to the layer's dtype of float32, which is new behavior in TensorFlow 2.  The layer has dtype float32 because its dtype defaults to floatx.

If you intended to run this layer in float32, you can safely ignore this warning. If in doubt, this warning is likely only an issue if you are porting a TensorFlow 1.X model to TensorFlow 2.

To change all layers to have dtype float64 by default, call `tf.keras.backend.set_floatx('float64')`. To change just this layer, pass dtype='float64' to the layer constructor. If you are the author of this layer, you can disable autocasting by passing autocast=False to the base Layer constructor.


Epoch: 0, loss:3.8782,  mae:1.3453,  val_loss:5.2470,  val_mae:1.4412,  
....................................................................................................
Epoch: 100, loss:0.0652,  mae:0.1886,  val_loss:0.0886,  val_mae:0.2114,  
.......................................[...]<jupyter_text>Determining error<jupyter_code>plotter = tfdocs.plots.HistoryPlotter(smoothing_std=2)
plotter.plot({'Early Stopping': early_history}, metric = "mae")
plt.ylabel('MAE [\u03BCm]')
loss, mae = model.evaluate(test_set_n[FEATURES], test_labels, verbose=2)
print("Testing set Mean Abs Error: {:5.2f} \u03BCm".format(mae))<jupyter_output>7/7 - 0s - loss: 0.0197 - mae: 0.0937
Testing set Mean Abs Error:  0.09 μm<jupyter_text>Prediction<jupyter_code>test_predictions = model.predict(test_set_n[FEATURES]).flatten() # using test set to make predictions
a = plt.axes(aspect='equal')
plt.scatter(test_labels, test_predictions)
plt.xlabel('True Values [\u03BCm]')
plt.ylabel('Predictions [\u03BCm]')
lims = [0, 1.5]
plt.xlim(lims)
plt.ylim(lims)
_ = plt.plot(lims, lims)
error =[a_i - b_i for a_i, b_i in zip(test_predictions, test_labels)]
plt.hist(error, bins = 20, range = [-0.4,0.4])
plt.xlabel("Prediction Error [\u03BCm]")
_ = plt.ylabel("Count")<jupyter_output><empty_output><jupyter_text>Eliminating obvious outliers in test set if any<jupyter_code># Simple anomaly detection function
# Takes a list of numbers and a threshold value. Removes any number that is a threshold times standard deviation away from the mean from the list
thresh = 3 # should be adjusted for different test data size
def delouts(mylist,thresh):
  num = 1
  while (num !=0 ):
    num = 0;
    std = np.std(mylist)
    mean = np.mean(mylist)
    for val in mylist:
      if (np.abs(mean-val)>thresh*std):
        mylist.remove(val)
        num = 1;
print("Error before removing outliers (linear model): "+str(round(np.std(error),3))+" [\u03BCm]")
delouts(error,thresh)
print("Error after removing outliers (if any): "+str(round(np.std(error),3))+" [\u03BCm]")
print("Error before correcting for size: "+str(round(np.std(test_labels),3))+" [\u03BCm]")<jupyter_output>Error before removing outliers (linear model): 0.14 [μm]
Error after removing outliers (if any): 0.104 [μm]
Error before correcting for size: 0.322 [μm]<jupyter_text>Training with neural network<jupyter_code>model = build_neural_network() 
EPOCHS = 2000
early_stop = keras.callbacks.EarlyStopping(monitor='val_loss', patience=10) # stop training whenever validation set stops showing improvement
early_history = model.fit(training_set_n[FEATURES], training_labels, 
                    epochs=EPOCHS, validation_split = 0.3, verbose=0, 
                    callbacks=[early_stop, tfdocs.modeling.EpochDots()])<jupyter_output>WARNING:tensorflow:Layer dense_17 is casting an input tensor from dtype float64 to the layer's dtype of float32, which is new behavior in TensorFlow 2.  The layer has dtype float32 because its dtype defaults to floatx.

If you intended to run this layer in float32, you can safely ignore this warning. If in doubt, this warning is likely only an issue if you are porting a TensorFlow 1.X model to TensorFlow 2.

To change all layers to have dtype float64 by default, call `tf.keras.backend.set_floatx('float64')`. To change just this layer, pass dtype='float64' to the layer constructor. If you are the author of this layer, you can disable autocasting by passing autocast=False to the base Layer constructor.


Epoch: 0, loss:2.1002,  mae:1.0485,  val_loss:2.4174,  val_mae:1.0366,  
....................................................................................................
Epoch: 100, loss:0.0176,  mae:0.0855,  val_loss:0.0198,  val_mae:0.0946,  
.......................................[...]<jupyter_text>Determining error<jupyter_code>plotter = tfdocs.plots.HistoryPlotter(smoothing_std=2)
plotter.plot({'Early Stopping': early_history}, metric = "mae")
plt.ylabel('MAE [\u03BCm]')
loss, mae = model.evaluate(test_set_n[FEATURES], test_labels, verbose=2)
print("Testing set Mean Abs Error: {:5.2f} \u03BCm".format(mae))<jupyter_output>7/7 - 0s - loss: 0.0124 - mae: 0.0800
Testing set Mean Abs Error:  0.08 μm<jupyter_text>Prediction<jupyter_code>test_predictions = model.predict(test_set_n[FEATURES]).flatten() # using test set to make predictions
a = plt.axes(aspect='equal')
plt.scatter(test_labels, test_predictions)
plt.xlabel('True Values [\u03BCm]')
plt.ylabel('Predictions [\u03BCm]')
lims = [0, 1.5]
plt.xlim(lims)
plt.ylim(lims)
_ = plt.plot(lims, lims)
error =[a_i - b_i for a_i, b_i in zip(test_predictions, test_labels)]
plt.hist(error, bins = 20, range = [-0.4,0.4])
plt.xlabel("Prediction Error [\u03BCm]")
_ = plt.ylabel("Count")<jupyter_output><empty_output><jupyter_text>Eliminating obvious outliers in test set if any<jupyter_code>print("Error before removing outliers (neural network): "+str(round(np.std(error),3))+" [\u03BCm]")
delouts(error,thresh)
print("Error after removing outliers (if any): "+str(round(np.std(error),3))+" [\u03BCm]")
print("Error before correcting for size: "+str(round(np.std(test_labels),3))+" [\u03BCm]")<jupyter_output>Error before removing outliers (neural network): 0.111 [μm]
Error after removing outliers (if any): 0.091 [μm]
Error before correcting for size: 0.322 [μm]
{'fraction_non_alphanumeric': 0.13511079084849575, 'fraction_numerical': 0.029634300126103404, 'mean_word_length': 4.673479816044967, 'pattern_counts': {'":': 0, '<': 61, '<?xml version=': 0, '>': 62, 'https://': 2, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '17', 'dir': 'jupyter-structured-clean-dedup', 'id': '489697', 'license': 'MIT', 'n_tokens_mistral': 3898, 'n_tokens_neox': 3556, 'n_words': 1050, 'path': 'Keras_Regression_on_beads.ipynb', 'repo_name': 'masterika/Fluorescence_BSVC', 'size': '11102'}