{ "cells": [ { "cell_type": "markdown", "metadata": { "_cell_guid": "d9837733-82a0-65f4-71a3-7e9ddbe669ae" }, "source": [ "This kernel is intended as a simple end-to-end tutorial with local validation using **mxnet**, one of the few deep learning frameworks for R. There are not that many mxnet examples in the web so I hope this template will be a usefull starting point for those beginning to build their first convolutional networks. Mxnet updated version has more functions to explore that I will not be using, exploring them can be your next step.\n", "\n", "Being this just a template-tutorial and limited by Kaggle kernel environment we will be using a very much reduced image resolution 64x64 and cpu mode (convnets training is much faster on gpus, more on that later). Also convnet architecture will be very simple, just two convolutional layers with pooling. \n", "\n", "\n", "With this in mind, lets make clear what the kernel will and achieve:\n", "\n", "- A simple image preprocessing pipeline that ends with dimensions ready for mxnet training.\n", "- A two convolutional layer, two pooling net for clasification.\n", "- A clasification result almost as good as one achieved by a **dart throwing monkey**, (with three classes a monkey would get a logloss of (log(3) = 1.098612))\n", "\n", "*( **EDIT**: We get a probably too optimistic error in our local validation of 0.9539. Better than expected. **A submission with this exact parameters gets 1.04438 on LB**. Worse than sample submission benchmark but still better than a monkey! )*\n", "\n", "Some parts of this code are taken from this fantastic kernel for State Farm Kaggle Competition [here](https://www.kaggle.com/gmilosev/r-mxnet) \n", "\n", "If you find this kernel useful remember to upvote! " ] }, { "cell_type": "code", "execution_count": null, "metadata": { "_cell_guid": "7928c643-fbcc-72b7-5ebc-daa15f3b1dd0" }, "outputs": [], "source": [ "library(dplyr)\n", "library(EBImage)\n", "library(mxnet)\n", "library(nnet)" ] }, { "cell_type": "markdown", "metadata": { "_cell_guid": "92c216a0-80f8-9db4-e1e6-35c378c1fdcb" }, "source": [ "Ok, let's preprocess images:" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "_cell_guid": "a0f2a9ba-9d2d-b207-2ecb-0d9985fc5528" }, "outputs": [], "source": [ "execute <- FALSE #KERNEL IS KILLED SO WILL NOT RUN IN KAGGLE ENVIRONMENT, SET TO TRUE LOCALLY\n", "\n", "if(execute == T){\n", " \n", "paths <- c(\"../input/train/Type_1\", \"../input/train/Type_2\", \"../input/train/Type_3\")\n", "\n", "\n", "\n", "in_type_counter <- 1 \n", "\n", "for (t in paths){\n", " \n", "patients <- dir(t)\n", "\n", "#For this simple example We rescale photos to 64*64 (= 4096)\n", "\n", "\n", "n_columnas <- 1 + 1 + 12288\n", "\n", "ordered_images <- data.frame(matrix(nrow = length(patients), ncol = n_columnas))\n", "colnames(ordered_images) <- c(\"paciente\", \"Type\", paste(\"R\", (1:4096), sep = \"\"), paste(\"G\", (1:4096), sep = \"\"), paste(\"B\", (1:4096), sep = \"\"))\n", "\n", "contador <- 1\n", "\n", "\n", "if(in_type_counter == 1){\n", "ordered_images$Type <- 1 #(*) ojo , asignaci\u00f3n manual segun el folder que est\u00e9 procesando\n", "}\n", "\n", "if(in_type_counter == 2){\n", "ordered_images$Type <- 2 #(*) ojo , asignaci\u00f3n manual segun el folder que est\u00e9 procesando\n", "}\n", "\n", "if(in_type_counter == 3){\n", "ordered_images$Type <- 3 #(*) ojo , asignaci\u00f3n manual segun el folder que est\u00e9 procesando\n", "}\n", "\n", "\n", "\n", "#+++\n", "mom_inicio <- Sys.time()\n", "print(\"beginning calculation: \")\n", "print(mom_inicio)\n", "#+++\n", "\n", "\n", "#THIS FOR LOOP WILL OPEN IMAGE, RESCALE IT TO 64X64 AND ASSIGN PIXELS INFO TO A DATA FRAME, ONE PATIEN PER ROW\n", "# NOTICE THAT WE ARE DISTORTING IMAGE PROPORTIONS RESIZING TO 64X64. IT IS A FAST WAY TO OBTAIN ALL IMAGES OF SAME SIZE, BETTER VERSIONS SHOULD CROP INSTEAD.\n", "\n", "\n", "#(YOU CAN ADJUST DIMENSIONS TO BIGGER SIZE WHEN RUNNING LOCALY)+++++++++++++++++++++++++++++++++++++++++++++++\n", "\n", "for (p in patients){\n", " \n", " cat(\"contador:\", contador, \" paciente:\", p, \"\\n\")\n", " \n", " \n", "\n", " ordered_images$paciente[contador] <- p\n", "\n", "\n", " imagen_paciente <- readImage(paste(t, p, sep = \"/\")) #abre imagen de cada paciente...\n", " imagen_paciente <- resize(imagen_paciente, w = 64, h = 64)# ((*)!!! OJO, NO CONSERVARA LA PROPORCION, se deforma, OPCION SOLO PARA DEMOSTRACION Y/O BASELINE)\n", "\n", "\n", " ordered_images[contador, c(3:4098)] <- imagen_paciente[, , 1]\n", " ordered_images[contador, c(4099:8194)] <- imagen_paciente[, , 2]\n", " ordered_images[contador, c(8195:12290)] <- imagen_paciente[, , 3]\n", " \n", "contador <- contador + 1 \n", " \n", "}\n", "\n", "#+++\n", "print(\"calculation time: \")\n", "print(Sys.time() - mom_inicio)\n", "#+++\n", "\n", "\n", "if(in_type_counter == 1){\n", "ordered_images_Type_1 <- ordered_images\n", "}\n", "\n", "if(in_type_counter == 2){\n", "ordered_images_Type_2 <- ordered_images\n", "}\n", "\n", "if(in_type_counter == 3){\n", "ordered_images_Type_3 <- ordered_images\n", "}\n", "\n", "\n", "in_type_counter <- in_type_counter + 1 \n", "}\n", "\n", "\n", "\n", "\n", "ordered_images_all <- bind_rows(ordered_images_Type_1, ordered_images_Type_2, ordered_images_Type_3)\n", " \n", " }" ] }, { "cell_type": "markdown", "metadata": { "_cell_guid": "a8c4470a-e951-17ca-2b70-f0953cb32e55" }, "source": [ "We now have **a single data frame containing one patient per row**. First two columns are patient id and cancer type, the rest are pixel intensities for the three colour channels, R G and B. \n", "Lets just check the first photo, just to make sure that image info was correctly stored:" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "_cell_guid": "b72631f7-2529-e919-2bd1-1b79e1f9ecff" }, "outputs": [], "source": [ "if(execute == T){\n", "\n", "row <- 1\n", "\n", "img_comprobacion <- Image(ordered_images_all[row, c(3:12286)], dim = c(64, 64, 3), colormode = 2)\n", "display(img_comprobacion, method = \"raster\")\n", "\n", "img_comprobacion\n", " \n", "}" ] }, { "cell_type": "markdown", "metadata": { "_cell_guid": "863d5e74-8273-395c-ffa6-a18f2248fa26" }, "source": [ "So that worked ok. Notice the low resolution and also the dimension distortion. We have made al images to be squares distorting image proportions what is not best practice. Better image processing should frequently relay more on cropping than on blind resizing. For our simple template though, it will be good enough.\n", "\n", "I have already warned that this convnet will not learn much by now. One reason is image resolution. Another reason of great importance is the amount of data we are feeding the net with. How many observations does a convnet need to learn? The answer is data specific and also depends on the architeture of the net. \n", "If I were asked how much data do I need for a \"not too complex\"\" problem with a convnet, (provided that I have the computing power), I would begin by asking 100.000 observations. That is my \"rule of thumb\" for a \"decent \" minimum amount of data for deep learning. In this case we obviously are not even close.\n", "\n", "What takes us to another important preprocessing step with convnets. **Image augmentation**. You can and you must either get more data or create it yourself modifying what you have. There you can use small scale changes, proportion distortion, rotations... We are not going to do anything of that here but keep in mind in this case it is a must given the sample size.\n", "\n", "And, connected to this... we have a quite unbalanced dataset. That is usually a problem for clasifications with many models and it is definitely a problem for convnets if we don't do some kind of rebalancing previously. (Undersampling the most frequent class, oversampling the others, etc.)\n", "\n", "For our template we are just going to undersample the most frequent two classes and get a ridiculously small, but balanced, dataset for our convnet:\n", "It is a good idea to plan previously to image augmentation how this class balance will be achieved so that image augmentation phase can obtain enough data of each class.\n", "\n", "Let's have a look at target distribution:" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "_cell_guid": "c168c678-2983-0bbb-ddd4-dfc0cc3ab1b7" }, "outputs": [], "source": [ "if(execute == T){\n", "#target distribution:\n", "table(ordered_images_all$Type)\n", " \n", "# **IMPORTANT STEP** lets transform label beginning with class 0, necessary for mxnet:\n", "ordered_images_all$Type <- ordered_images_all$Type - 1\n", "\n", " \n", "}" ] }, { "cell_type": "markdown", "metadata": { "_cell_guid": "1bb1aa73-7d3f-3226-ef67-b5778540c77e" }, "source": [ "First step: we sample 75 observations, 25 of each class, for our validation set:" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "_cell_guid": "1efacb5d-61e0-b482-8adf-08bb3248021f" }, "outputs": [], "source": [ "if(execute == T){\n", "\n", "set.seed(9)\n", "set_validacion <- ordered_images_all %>% \n", " group_by(Type) %>% \n", " sample_n(25) %>% \n", " ungroup()\n", " \n", "}" ] }, { "cell_type": "markdown", "metadata": { "_cell_guid": "b74e9f66-5767-49fb-37a4-89cb31c87b34" }, "source": [ "Now we undersample types 2 and 3, and use all type 1 to get a balanced train set of 3*225 cases. And we shuffle that train set (We want our min batches with a distribution as close as possible to train)" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "_cell_guid": "95a19f1b-5b87-e96c-86d8-3445edb2efb5" }, "outputs": [], "source": [ "if(execute == T){\n", "\n", "set_train_unbalanced <- ordered_images_all[ordered_images_all$paciente %in% setdiff(ordered_images_all$paciente, set_validacion$paciente), ]\n", "\n", "#primera muestra train undersampled:\n", "set.seed(9)\n", "undersample_1_set_train <- set_train_unbalanced %>% \n", " group_by(Type) %>% \n", " sample_n(225) %>% \n", " ungroup() %>% \n", " sample_frac(1) %>% \n", " sample_frac(1) #WE SHUFFLE TWICE\n", "\n", "# dim(undersample_1_set_train)\n", "# #[1] 675 12290\n", "# dim(set_validacion)\n", "# #[1] 75 12288\n", "\n", "\n", "#_____train and validaci\u00f3n labels___________________________\n", "target_undersample_1 <- undersample_1_set_train$Type\n", "label_set_validacion <- set_validacion$Type\n", "#____________________________________________________________\n", "\n", "\n", "#_____train set1(undersampled) y and validation set, both balanced_______\n", "undersample_1_set_train <-undersample_1_set_train[, c(3:12290)]\n", "set_validacion <-set_validacion[, c(3:12290)]\n", "#_____target y label de validaci\u00f3n_________________________________________\n", "\n", "\n", "# AND NOW WE HAVE TO BE CAREFUL WITH DIMENSIONS, FIRST TRASPOSE, THEN DIMENSIONS:\n", "undersample_1_set_train <- t(undersample_1_set_train)\n", "dim(undersample_1_set_train) <- c(64, 64, 3, 675)\n", "\n", "\n", "set_validacion <- t(set_validacion)\n", "dim(set_validacion) <- c(64, 64, 3, 75)\n", " \n", " }" ] }, { "cell_type": "markdown", "metadata": { "_cell_guid": "8c12c204-c648-1b2a-10ed-eb704f137eab" }, "source": [ "To the convnet now. (I am assuming some basic knowledge on convolutional nets here)" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "_cell_guid": "886eafdb-e157-fd06-ee64-31ac4f33fea9" }, "outputs": [], "source": [ "if(execute == T){\n", "\n", "#LETS CHOOSE SOME PARAMETTERS IN OUR CONVNET\n", "#Just some parametters that perform not too badly:\n", "n_output <- 3 # we have three classes, will use it for calculating number of nodes in fc1\n", "num_filters_conv2 = 14\n", "num.round = 101\n", "learning.rate = 0.1\t\n", "momentum = 0.0056\n", "weight_decay = 0.0046\t\n", "initializer = 0.0667\t\n", "\n", "\n", "# AND SOME HELPER FUNCTIONS:\n", "mLogLoss.normalize = function(p, min_eta=1e-15, max_eta = 1.0){\n", " #min_eta\n", " for(ix in 1:dim(p)[2]) {\n", " p[,ix] = ifelse(p[,ix]<=min_eta,min_eta,p[,ix]);\n", " p[,ix] = ifelse(p[,ix]>=max_eta,max_eta,p[,ix]);\n", " }\n", " #normalize\n", " for(ix in 1:dim(p)[1]) {\n", " p[ix,] = p[ix,] / sum(p[ix,]);\n", " }\n", " return(p);\n", "}\n", "\n", "# helper function\n", "#calculates logloss\n", "mlogloss = function(y, p, min_eta=1e-15,max_eta = 1.0){\n", " class_loss = c(dim(p)[2]);\n", " loss = 0;\n", " p = mLogLoss.normalize(p,min_eta, max_eta);\n", " for(ix in 1:dim(y)[2]) {\n", " p[,ix] = ifelse(p[,ix]>1,1,p[,ix]);\n", " class_loss[ix] = sum(y[,ix]*log(p[,ix]));\n", " loss = loss + class_loss[ix];\n", " }\n", " #return loss\n", " return (list(\"loss\"=-1*loss/dim(p)[1],\"class_loss\"=class_loss));\n", "}\n", "\n", "# mxnet specific logloss metric\n", "mx.metric.mlogloss <- mx.metric.custom(\"mlogloss\", function(label, pred){\n", " p = t(pred);\n", " m = mlogloss(class.ind(label),p);\n", " gc();\n", " return(m$loss);\n", "})\n", " \n", "}" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "_cell_guid": "b3238f55-b6ab-19bb-e8fc-116ecd752736" }, "outputs": [], "source": [ "if(execute == T){\n", "\n", "train.x <- undersample_1_set_train\n", "train.y <-target_undersample_1\n", "\n", "\n", "#____________\n", "data = mx.symbol.Variable('data')\n", "#FIRST CONVOLUITIONAL LAYER + POOLING\n", "conv1 = mx.symbol.Convolution(data=data, kernel=c(3, 3), num_filter = 3) \n", "relu1 = mx.symbol.Activation(data=conv1, act_type=\"relu\") \n", "pool1 = mx.symbol.Pooling(data=relu1, pool_type=\"max\", kernel=c(2,2), stride=c(2,2))\n", "\n", "\n", "#SECOND CONVOLUTIONAL LAYER + POOLING\n", "conv2 = mx.symbol.Convolution(data=pool1, kernel=c(3,3), num_filter = num_filters_conv2) \n", "relu2 = mx.symbol.Activation(data=conv2, act_type=\"relu\")\n", "pool2 = mx.symbol.Pooling(data=relu2, pool_type=\"max\",kernel=c(2,2), stride=c(2,2)) \n", "\n", "#FLATTEN THE OUTPUT\n", "flatten = mx.symbol.Flatten(data=pool2) \n", "\n", "#FEED FULLY CONNECTED LAYER, NUMBER OF HIDDEN NODES JUST GEOMETRIC MEAN OF INPUT(14.5 * 14.5 * 13 = 2733.25) AND OUTPUT (3), sqrt(2733.25*3) = 91\n", "input_previo_a_filtroconv2 <- 14.5*14.5\n", "n_input <- input_previo_a_filtroconv2 * num_filters_conv2\n", "num_hidden_fc1 <- round(sqrt(n_input*n_output)) \n", "\n", "fc1 = mx.symbol.FullyConnected(data=flatten, num_hidden=84) \n", "relu4 = mx.symbol.Activation(data=fc1, act_type=\"relu\") \n", "\n", "#____________\n", "\n", "\n", "fc2 = mx.symbol.FullyConnected(data=relu4, num_hidden=3) #ESTA PARA CLASIFICACION\n", "\n", "mi_softmax = mx.symbol.SoftmaxOutput(data=fc2)\n", "\n", "\n", "devices <- mx.cpu()\n", "#(IF MXNET IS COMPILED TO WORK WITH GPU YOU CHANGE THE FOLLOWING LINE FOR MUCH HIGHER SPEED)\n", "#devices <- mx.gpu()\n", "\n", "\n", "\n", "mx.set.seed(0)\n", "\n", "#___\n", "tic <- proc.time()\n", "#___\n", "model <- mx.model.FeedForward.create( mi_softmax #for clasification\n", " , X=train.x\n", " , y=train.y\n", " , eval.data = list(\"data\" = set_validacion,\"label\" = label_set_validacion) \n", " , ctx=devices\n", " , num.round=num.round\n", " , array.batch.size = 75 \n", " , learning.rate = learning.rate\n", " , momentum = momentum\n", " , wd=weight_decay\n", " , eval.metric = mx.metric.mlogloss \n", " , initializer=mx.init.uniform(initializer) \n", " #, epoch.end.callback = mx.callback.save.checkpoint(\"modelo_guardado_ccs\") #(TO SAVE MODEL AT EVERY ITERATION)\n", " , batch.end.callback = mx.callback.log.train.metric(10)#, log) \n", " , array.layout=\"columnmajor\"\n", " ) \n", "\n", "#___\n", "print(proc.time() - tic)\n", "#___\n", " \n", " } \n", "\n", "# [101] Train-mlogloss=0.904353529359042\n", "# [101] Validation-mlogloss=0.953891121726141\n", "\n" ] }, { "cell_type": "markdown", "metadata": { "_cell_guid": "d3d0d818-3e37-f151-35f0-9a671ac2b36d" }, "source": [ "We get a probably too optimistic error in our local validation of 0.9539. Better than expected. (A submission with this exact parameters gets 1.04438 on LB. worse than sample submission benchmark but still better than a monkey!). \n", "Also notice we have not used the classes balance in our training so the performance of the model has a bit more merit than it seems.\n", "\n", "If we want we can predict on our validation set, and have a look at probabilities predicted vs. correct label:" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "_cell_guid": "266ec0c2-8ba6-309b-adc1-c828c076dbe6" }, "outputs": [], "source": [ "if(execute == T){\n", "\n", "preds <- predict(model, set_validacion,\n", " ctx = NULL,\n", " array.layout = \"auto\")\n", "\n", "#WE CAN INSPECT OUR LIMITED VALIDATION SET, PROBABILITIES:\n", "predicciones <- t(preds)\n", "predicciones <- cbind(label_set_validacion, predicciones)\n", "head(predicciones)\n", " \n", "}" ] }, { "cell_type": "markdown", "metadata": { "_cell_guid": "068d04e3-fe8f-eb69-af6d-aeba904bba21" }, "source": [ "Ok, so this is a template. Once you understand how this very basic example works I suggest a list of possible steps:\n", "\n", " * Add more convolutional layers\n", " * Try different number of filters on each layer as you advance through layers, keep in mind the output size you are building (Pooling reduces this size, filters increase it)\n", " * Keep in mind also that your output will affect the optimal number of nodes in the fully conected layer.\n", " * You can use dropping layers just like in any other NN.\n", " * Add or remove pooling layers depending on how you want your kernels learn across the scale of your image.\n", " * Try different kernel sizes.\n", " * Explore how minibatch size affects training.\n", " \n", "After all this you can make a more realistic approach trying bigger images and MANY more images through augmentation so that the net can really learn.\n", "For this you will have to make mxnet work in GPU mode so that training takes much less time. \n", "\n", "So, I hope this was of some use, **if you liked it remember you can upvote in the upper part of the page**, what is always motivating! :-)" ] } ], "metadata": { "_change_revision": 0, "_is_fork": false, "kernelspec": { "display_name": "R", "language": "R", "name": "ir" }, "language_info": { "codemirror_mode": "r", "file_extension": ".r", "mimetype": "text/x-r-source", "name": "R", "pygments_lexer": "r", "version": "3.4.0" } }, "nbformat": 4, "nbformat_minor": 0 }