{"cells":[{"metadata":{"_uuid":"635b9590de5790f9a4f6a5aec819f834e8bf5fb6"},"cell_type":"markdown","source":"# KERAS - TRAINING WITH FLOAT16 - Test kernel 2\n\n## Introduction\n\nDue to the big size of the images and the required detail for a good model to work in the Human Protein Atlas Image Classification Challenge, one of the possibilities of reducing the amount of memory needed is training with `float16` precision.    \n\nIn Keras, this should be done simply by setting `K.set_floatx('float16')`, however a few other things must be done to avoid `nan` values and to use the BatchNormalization layer using Tensorflow backend, which requires `float32` in all cases. (This Kernel was not tested for other backends and they may work differently)   \n\nAlso, after fixing the normalization layer, it will be necessary to fix the optimizer for conflicting types. \n\n### Warning:\n\nAlthough this kernel shows that it's possible to find a batch size that works in `float16` while not in `float32`, it seems slower and occupies more memory than [Test Kernel 1](https://www.kaggle.com/danmoller/keras-training-with-float16-test-kernel-1), which on its side doesn't seem to favor `float16`. \n\nBut these tests are very fresh and if you find a bug or can solve the issue, please let us know :)   \n\n\n## Differences between Test Kernel 1 and Test Kernel 2\n\nIn this kernel 2, we change the batch normalization layer to use `float16`. \n  \n - This skips using Tensorflow's fused batch normalization and uses a regular batch normalization    \n - This also skips the need of changing the optimizer for different formats   \n \n In [Test Kernel 1](https://www.kaggle.com/danmoller/keras-training-with-float16-test-kernel-1), we use the original batch normalization operations from Keras (they are faster and occupy less memory)    \n\n\n## Getting started - Setting to float16 and avoiding NaNs\n\nIn order to do this, simply run these commands before anything else.    \nThe purpose of setting `epsilon` to a bigger value is because the default value is too little for `float16` and will cause `nan` loss values during training. \n\n**The epsilon value is not optimized, you may want to test it with 1e-3 or other values**"},{"metadata":{"trusted":true,"_uuid":"f7acd238c18d631aed8b5d38e3e47fb4bedef540"},"cell_type":"code","source":"import keras.backend as K\nK.set_floatx('float16')\nK.set_epsilon(1e-4) #default is 1e-7\n\nisTestMode = False #use this for training short epochs to quickly see the results\nepochs = 2 \nbatchSize = 168","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"59d060743cc3b8639856512fd0b7c15991be7cb4"},"cell_type":"markdown","source":"## Fixing the BatchNormalization layer\n\nBecause of Tensorflow's requirement of using `float32` in batch normalization, the setting above will break some things because Keras will send `float16` values to Tensorflow.\n\nThus, we will create custom weight initializers and a custom BatchNormalization layer:"},{"metadata":{"trusted":true,"_uuid":"c7e66838182ac2c9302e8e49fdfa1c50bc116091"},"cell_type":"code","source":"from keras.layers import BatchNormalization, Layer\nfrom keras.initializers import Initializer\nfrom keras.backend.tensorflow_backend import tf, _regular_normalize_batch_in_training\n\n\n#custom initializers to force float32\nclass Ones32(Initializer):\n    def __call__(self, shape, dtype=None):\n        return K.constant(1, shape=shape, dtype='float32')\n\nclass Zeros32(Initializer):\n    def __call__(self, shape, dtype=None):\n        return K.constant(0, shape=shape, dtype='float32')\n    \n\n\nclass BatchNormalizationF16(Layer):\n\n    def __init__(self,\n                 axis=-1,\n                 momentum=0.99,\n                 epsilon=1e-3,\n                 center=True,\n                 scale=True,\n                 beta_initializer='zeros',\n                 gamma_initializer='ones',\n                 moving_mean_initializer='zeros',\n                 moving_variance_initializer='ones',\n                 beta_regularizer=None,\n                 gamma_regularizer=None,\n                 beta_constraint=None,\n                 gamma_constraint=None,\n                 **kwargs):\n        super(BatchNormalizationF16, self).__init__(**kwargs)\n        self.supports_masking = True\n        self.axis = axis\n        self.momentum = momentum\n        self.epsilon = epsilon\n        self.center = center\n        self.scale = scale\n        self.beta_initializer = initializers.get(beta_initializer)\n        self.gamma_initializer = initializers.get(gamma_initializer)\n        self.moving_mean_initializer = initializers.get(moving_mean_initializer)\n        self.moving_variance_initializer = (\n            initializers.get(moving_variance_initializer))\n        self.beta_regularizer = regularizers.get(beta_regularizer)\n        self.gamma_regularizer = regularizers.get(gamma_regularizer)\n        self.beta_constraint = constraints.get(beta_constraint)\n        self.gamma_constraint = constraints.get(gamma_constraint)\n\n    def build(self, input_shape):\n        dim = input_shape[self.axis]\n        if dim is None:\n            raise ValueError('Axis ' + str(self.axis) + ' of '\n                             'input tensor should have a defined dimension '\n                             'but the layer received an input with shape ' +\n                             str(input_shape) + '.')\n        self.input_spec = InputSpec(ndim=len(input_shape),\n                                    axes={self.axis: dim})\n        shape = (dim,)\n\n        if self.scale:\n            self.gamma = self.add_weight(shape=shape,\n                                         name='gamma',\n                                         initializer=self.gamma_initializer,\n                                         regularizer=self.gamma_regularizer,\n                                         constraint=self.gamma_constraint)\n        else:\n            self.gamma = None\n        if self.center:\n            self.beta = self.add_weight(shape=shape,\n                                        name='beta',\n                                        initializer=self.beta_initializer,\n                                        regularizer=self.beta_regularizer,\n                                        constraint=self.beta_constraint)\n        else:\n            self.beta = None\n        self.moving_mean = self.add_weight(\n            shape=shape,\n            name='moving_mean',\n            initializer=self.moving_mean_initializer,\n            trainable=False)\n        self.moving_variance = self.add_weight(\n            shape=shape,\n            name='moving_variance',\n            initializer=self.moving_variance_initializer,\n            trainable=False)\n        self.built = True\n\n    def call(self, inputs, training=None):\n        input_shape = K.int_shape(inputs)\n        # Prepare broadcasting shape.\n        ndim = len(input_shape)\n        reduction_axes = list(range(len(input_shape)))\n        del reduction_axes[self.axis]\n        broadcast_shape = [1] * len(input_shape)\n        broadcast_shape[self.axis] = input_shape[self.axis]\n\n        # Determines whether broadcasting is needed.\n        needs_broadcasting = (sorted(reduction_axes) != list(range(ndim))[:-1])\n\n        def normalize_inference():\n            if needs_broadcasting:\n                # In this case we must explicitly broadcast all parameters.\n                broadcast_moving_mean = K.reshape(self.moving_mean,\n                                                  broadcast_shape)\n                broadcast_moving_variance = K.reshape(self.moving_variance,\n                                                      broadcast_shape)\n                if self.center:\n                    broadcast_beta = K.reshape(self.beta, broadcast_shape)\n                else:\n                    broadcast_beta = None\n                if self.scale:\n                    broadcast_gamma = K.reshape(self.gamma,\n                                                broadcast_shape)\n                else:\n                    broadcast_gamma = None\n                return tf.nn.batch_normalization(#K.batch_normalization(\n                    inputs,\n                    broadcast_moving_mean,\n                    broadcast_moving_variance,\n                    broadcast_beta,\n                    broadcast_gamma,\n                    #axis=self.axis,\n                    self.epsilon)#epsilon=self.epsilon)\n            else:\n                return tf.nn.batch_normalization(#K.batch_normalization(\n                    inputs,\n                    self.moving_mean,\n                    self.moving_variance,\n                    self.beta,\n                    self.gamma,\n                    #axis=self.axis,\n                    self.epsilon)#epsilon=self.epsilon)\n\n        # If the learning phase is *static* and set to inference:\n        if training in {0, False}:\n            return normalize_inference()\n\n        # If the learning is either dynamic, or set to training:\n        normed_training, mean, variance = _regular_normalize_batch_in_training(#K.normalize_batch_in_training(\n            inputs, self.gamma, self.beta, reduction_axes,\n            epsilon=self.epsilon)\n\n        if K.backend() != 'cntk':\n            sample_size = K.prod([K.shape(inputs)[axis]\n                                  for axis in reduction_axes])\n            sample_size = K.cast(sample_size, dtype=K.dtype(inputs))\n\n            # sample variance - unbiased estimator of population variance\n            variance *= sample_size / (sample_size - (1.0 + self.epsilon))\n\n        self.add_update([K.moving_average_update(self.moving_mean,\n                                                 mean,\n                                                 self.momentum),\n                         K.moving_average_update(self.moving_variance,\n                                                 variance,\n                                                 self.momentum)],\n                        inputs)\n\n        # Pick the normalized form corresponding to the training phase.\n        return K.in_train_phase(normed_training,\n                                normalize_inference,\n                                training=training)\n\n    def get_config(self):\n        config = {\n            'axis': self.axis,\n            'momentum': self.momentum,\n            'epsilon': self.epsilon,\n            'center': self.center,\n            'scale': self.scale,\n            'beta_initializer': initializers.serialize(self.beta_initializer),\n            'gamma_initializer': initializers.serialize(self.gamma_initializer),\n            'moving_mean_initializer':\n                initializers.serialize(self.moving_mean_initializer),\n            'moving_variance_initializer':\n                initializers.serialize(self.moving_variance_initializer),\n            'beta_regularizer': regularizers.serialize(self.beta_regularizer),\n            'gamma_regularizer': regularizers.serialize(self.gamma_regularizer),\n            'beta_constraint': constraints.serialize(self.beta_constraint),\n            'gamma_constraint': constraints.serialize(self.gamma_constraint)\n        }\n        base_config = super(BatchNormalizationF16, self).get_config()\n        return dict(list(base_config.items()) + list(config.items()))\n\n    def compute_output_shape(self, input_shape):\n        return input_shape\n\n","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"d5b2d8dd350d0ec02891b57b8b8bbf82cdfe0480"},"cell_type":"markdown","source":"## Fixing the optimizers\n\nIn this test kernel 2, we won't need to fix the optimizers, since everything is `float16`, there will be no mixing in the optimizer and Keras will handle everything properly."},{"metadata":{"trusted":true,"_uuid":"fa6f4a2c6e871238ab800682d0baa65f524e5bd3"},"cell_type":"code","source":"from keras.optimizers import SGD","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"a0ea277d245833f2593edb79b6a279372a55ef42"},"cell_type":"markdown","source":"# Training\n\nNow we're testing our changes in a simple model (this model is not really well thought for this competition, just an example). \n\nLets create the model after a few definitions for loading data and organizing it in images with 4 channels (RGBY). Don't forget to convert your data to `float16`."},{"metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true},"cell_type":"code","source":"competitionFolder = '../input/' #human-protein-atlas-image-classification/'\ntrainFolder = competitionFolder + 'train/'\ntestFolder = competitionFolder + 'test/'\nnClasses = 28\nside=512\noriginalSide = 512\ncropSide = 256\n","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"36fdd5e72da513dc686f714ac85e15f8d7c65170"},"cell_type":"code","source":"import numpy as np\nfrom PIL import Image\nimport random\n\n%matplotlib inline\nimport matplotlib.pyplot as plt\n\ndef loadClasses():\n    trainFile = competitionFolder + 'train.csv'\n    filesAndClasses = list()\n    \n    with open(trainFile, 'r') as f:\n        _ = next(f)\n        for row in f:\n            fields = row.split(',')\n            file = fields[0]\n            \n            classesNp = np.zeros((nClasses,), dtype='float16')\n            classes = fields[1].split(' ')\n            for c in classes: classesNp[int(c)] = 1\n                \n            filesAndClasses.append((trainFolder + file, classesNp))\n    return filesAndClasses\n\ndef loadImage(file):\n    colors = ['_red.png', '_green.png', '_blue.png', '_yellow.png']\n    images = [Image.open(file + color) for color in colors]\n    return np.stack(images, axis=-1)\n\n#flips a batch of images, flipMode is an integer in range(8)\ndef flip(x, flipMode):\n    if flipMode in [4,5,6,7]:\n        x = np.swapaxes(x,1,2)\n    if flipMode in [1,3,5,7]:\n        x = np.flip(x,1)\n    if flipMode in [2,3,6,7]:\n        x = np.flip(x,2)\n        \n    return x\n\ndef inspect(x, name):\n    print(name + \": \", 'shape:', x.shape, 'min:', x.min(), 'max:',x.max())\n    \ndef plotChannels(img, minVal = 0, maxVal = 255):\n    fig, ax = plt.subplots(1, img.shape[-1], figsize=(20,10))\n    for i in range(img.shape[-1]):\n        ax[i].imshow(img[:,:,i], vmin = minVal, vmax= maxVal)\n        \n    plt.show()\n    \ndef competitionMetric(true,pred):\n    pred = K.cast(K.greater(pred,0.5), K.floatx())\n    \n    groundPositives = K.sum(true, axis=0) + K.epsilon()\n    correctPositives = K.sum(true * pred, axis=0)\n    predictedPositives = K.sum(pred, axis=0) + K.epsilon()\n\n    precision = correctPositives / predictedPositives\n    recall = correctPositives / groundPositives\n\n    m = (2 * precision * recall) / (precision + recall + K.epsilon())\n\n    return K.mean(m)","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"a19171e393f4862421a74fbaae5b6efe1f589617"},"cell_type":"markdown","source":"## Data generator with cropping and flipping\n\nIn order to train even faster, we're creating a data generator that crops from the 512x512 images.   \nIt was not studied if this crop may hide areas that contain the target proteins, but it's very probable that the protein be present in this crop if it's big enough.\n\nHere we will be training with crops of size 256x256\n\n"},{"metadata":{"trusted":true,"_uuid":"af11d46679c2b5a9a30baaa2285e3e425cb553d2"},"cell_type":"code","source":"from keras.utils import Sequence\nfrom random import shuffle\n\n#works with channels last\nclass ImageLoader(Sequence):\n    \n    #class creator, use generationMode = 'predict' for returning only images without labels\n        #when using 'predict', pass only a list of files, not files and classes\n    def __init__(self, filesAndClasses, batchSize, generationMode = 'train'):\n        \n        self.filesAndClasses = filesAndClasses\n        self.batchSize = batchSize\n        self.generationMode = generationMode\n        \n        assert generationMode in ['train', 'predict']\n            \n\n    #gets the number of batches this generator returns\n    def __len__(self):\n        l,rem = divmod(len(self.filesAndClasses), self.batchSize)\n        return (l + (1 if rem > 0 else 0))\n    \n    #shuffles data on epoch end\n    def on_epoch_end(self):\n        if self.generationMode == 'train':\n            shuffle(self.filesAndClasses)\n        \n    #gets a batch with index = i\n    def __getitem__(self, i):\n        \n        #x are images   \n        #y are labels\n        \n        pairs = self.filesAndClasses[i*self.batchSize:(i+1)*self.batchSize]\n        if self.generationMode == 'train':\n            files, classes = zip(*pairs) \n            y = np.stack(classes, axis=0)\n            x = [loadImage(f) for f in files]\n        elif self.generationMode == 'predict':\n            files = pairs\n            x = [loadImage(f) for f in files]\n        else:\n            raise Exception(\"ImageLoader does not support 'generationMode' of type \" + self.generationMode)\n    \n        x = np.stack(x, axis=0)\n        \n        #cropping and flipping when training\n        if self.generationMode == 'train':\n            \n            startH = random.randint(0,side - cropSide)\n            startW = random.randint(0,side - cropSide)\n            \n            x = x[:, startW:startW + cropSide, startH:startH + cropSide]\n            \n            flipMode = random.randint(0,7) #see flip functoin defined above\n            x = flip(x, flipMode)\n\n        if self.generationMode == 'predict':\n            return x\n        else:\n            return x, y\n        \n","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"27ec522b5828336af58658ecdad024f491678f08"},"cell_type":"markdown","source":"## Loading data and creating generator\n\nLet's load the data an make a quick inspection of the generator. "},{"metadata":{"trusted":true,"_uuid":"a04a217cbecfe2ce4a4bee9e06151f4a41e49de6","scrolled":false},"cell_type":"code","source":"valBatchSize = batchSize // 4\n\ntrainFiles = loadClasses()\n\n#creating a fold for validation\nfold = 0\ntestLen = len(trainFiles)//5\ntestStart = fold * testLen\ntestEnd = testStart + testLen\n\nvalFiles = trainFiles[testStart:testEnd]\ntrainFiles = trainFiles[0:testStart] + trainFiles[testEnd:]\n\nif isTestMode:\n    valFiles = valFiles[:5*valBatchSize]\n    trainFiles = trainFiles[:5*batchSize]\n\n#creating train and val generators\ntrainGenerator = ImageLoader(trainFiles, batchSize)\nvalGenerator = ImageLoader(valFiles, valBatchSize)\n\n#quick check\ntrainFileList, trainLabels = zip(*trainFiles)\npredictGenerator = ImageLoader(trainFileList, batchSize, generationMode='predict')\nfor i in range(5):\n    originalX = predictGenerator[i]\n    x,y = trainGenerator[i]\n    inspect(x,'images')\n    inspect(y,'labels')\n    print('unique y: ', np.unique(y))\n    print('original')\n    plotChannels(originalX[0])\n    print('cropped')\n    plotChannels(x[0])","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"37c1ed433f161aa2ac7e51adeb765375e91aefbe"},"cell_type":"markdown","source":"## Model Creator\n\nHere, a simple model using the custom layer and optimizer, similar to a ResNet. \n\n"},{"metadata":{"trusted":true,"_uuid":"e9892f392bb48b5f848c4b0abe229666c88e69db"},"cell_type":"code","source":"from keras.layers import *\nfrom keras.models import Model\n\ndef modelCreator(convFilters, denseFilters):\n                \n    ######################################## definitions for layers ###########################3\n    \n    def denseBN(inputs, filters, activation, name):\n        out = Dense(filters, name = name, use_bias=False)(inputs)\n        out = BatchNormalizationF16(name = name + \"BN\")(out)\n        out = Activation(activation, name = name + \"ACT\")(out)\n        \n        return out\n    \n    def convBN(inputs, filters, kernelSize, activation, name):\n        out = Conv2D(filters, kernelSize, name=name, padding='same', use_bias=False)(inputs)\n        out = BatchNormalizationF16(name = name + \"BN\")(out)\n        out = Activation(activation, name = name + 'ACT')(out)\n            \n        return out\n    \n    ##################################### block definitions ####################################\n\n    def downBlock(i, filters, inputs):\n        \n        name = str(i)\n        out = inputs\n        \n        #make maxpooling and resnet connection if not first block\n        if (i != 0):\n            \n            out = MaxPooling2D(poolSizes[i], name='Down' + name)(out)\n            connection = convBN(out, filters, 3, activation = 'linear', name = 'resConnDownB' + name)\n\n\n        out = convBN(out,filters, 3, activation='relu', name = 'downConvA' + name)\n        out = convBN(out, filters, 3, activation='relu', name = 'downConvB' + name )\n        out = convBN(out, filters, 3, activation='relu', name = 'downConvC' + name )\n        \n        #resnet connection\n        if i != 0:\n            out = convBN(out, filters, 3, activation = 'linear', name = 'resConnDownA' + name)\n            out = Add(name='resAddDown' + name)([out,connection])\n            out = BatchNormalizationF16(name = 'resNormDown' + name)(out)\n            out = Activation('relu', name = 'downAct' + name)(out)\n        \n        return out\n    \n    ####################################### model creation #############################################\n    \n        \n    poolSizes = [0,4,4,4,4]\n    \n    #notice we are training with 256x256 and validating with 512x512, thus None as size\n    inp = Input((None,None,4))\n    out = BatchNormalizationF16(name='initNorm')(inp)\n\n    for i,filts in enumerate(convFilters):\n        out = downBlock(i,filts,out)\n    \n    out = GlobalMaxPooling2D(name='globalPool')(out)\n    \n    for i, filts in enumerate(denseFilters):\n        out = denseBN(out, filts, 'relu', name = 'dense' + str(i))\n    \n    out = denseBN(out, nClasses, activation='sigmoid', name=\"FinalDense\")\n\n    model = Model(inp,out)\n                       \n    return model\n\n","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"cd9ba8458e15c1275eab726170c0e9d19f83840e"},"cell_type":"markdown","source":"## Creating, compiling and fitting\n\nLet's do it. \n\n**Warning:** the loss function selected may not be the best for this competition. "},{"metadata":{"trusted":true,"_uuid":"7e8ab71ef5b1a0dd49a184f1ce9c1607a896f99b"},"cell_type":"code","source":"model = modelCreator(convFilters =  [20,40,90,130,200],\n                     denseFilters = [100,50,30])\n\n#confirm dtype is float16\nprint(\"type is: \", K.dtype(model.get_layer('downConvA0').kernel))\n\n#use a regular SGD\nmodel.compile(optimizer = SGD(lr=0.01,momentum=.9), loss = 'categorical_crossentropy', metrics=[competitionMetric])\n\nmodel.fit_generator(trainGenerator,len(trainGenerator), \n                    validation_data = valGenerator, validation_steps = len(valGenerator),\n                   epochs = epochs, workers=5, max_queue_size=10)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"9e51189adc64e78baeb7bc468bed5fdfe1bf30a5"},"cell_type":"markdown","source":"# Saving and Loading\n\nIn order to save and load the model using custom layers and optimizers, one needs to create a custom objects dictionary to tell Keras how to recreate these objects.\n\nSo we should include our layer, initializers, custom metric and optimizer in this object."},{"metadata":{"trusted":true,"_uuid":"657dd3bf9a74d060583fc4263ecf552c6c92b8da"},"cell_type":"code","source":"from keras.models import load_model\n\ncustomObjects = {\n    'BatchNormalizationF16': BatchNormalizationF16,\n    'competitionMetric': competitionMetric,\n    'Ones32': Ones32,\n    'Zeros32': Zeros32\n}\n\nmodel.save('savedModel')\nloadedModel = load_model('savedModel', customObjects)\n\n#training starts from where it ended, including otimizer state\nloadedModel.fit_generator(trainGenerator,len(trainGenerator), \n                    validation_data = valGenerator, validation_steps = len(valGenerator),\n                   epochs = epochs, workers=5, max_queue_size=10)","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"6fa2a2678bf67bc66488b8cfc4d5a39079d0d552"},"cell_type":"markdown","source":"# Comparing with float32\n\nLet's try to train the same generator with same batch size on a model with precision `float32` and see the GPU return an \"out of memory (OOM)\" error.\n\nEven though, take a look at the warning at the beginning of this kernel and compare it with [Test Kernel 1](https://www.kaggle.com/danmoller/keras-training-with-float16-test-kernel-1)"},{"metadata":{"trusted":true,"_uuid":"e59f707486b19ca3e2c26cd95a47f798d0d99205"},"cell_type":"code","source":"dtype='float32'\nK.set_floatx(dtype)\n\n\nmodel = modelCreator(convFilters =  [20,40,90,130,200],\n                     denseFilters = [100,50,30])\n\n#confirm dtype is float32\nprint(\"type is: \", K.dtype(model.get_layer(\"downConvA0\").kernel))\n\n#use a regular SGD\nmodel.compile(optimizer = SGD(lr=0.01,momentum=.9), loss = 'categorical_crossentropy', metrics=[competitionMetric])\n\nmodel.fit_generator(trainGenerator,len(trainGenerator), \n                    validation_data = valGenerator, validation_steps = len(valGenerator),\n                   epochs = epochs, workers=5, max_queue_size=10)","execution_count":null,"outputs":[]}],"metadata":{"kernelspec":{"display_name":"Python 3","language":"python","name":"python3"},"language_info":{"name":"python","version":"3.6.6","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"}},"nbformat":4,"nbformat_minor":1}