ML4A_TSNE.ipynb

{
  "nbformat": 4,
  "nbformat_minor": 0,
  "metadata": {
    "kernelspec": {
      "display_name": "Python 3",
      "language": "python",
      "name": "python3"
    },
    "language_info": {
      "codemirror_mode": {
        "name": "ipython",
        "version": 3
      },
      "file_extension": ".py",
      "mimetype": "text/x-python",
      "name": "python",
      "nbconvert_exporter": "python",
      "pygments_lexer": "ipython3",
      "version": "3.5.2"
    },
    "colab": {
      "name": "ML4A_TSNE.ipynb",
      "private_outputs": true,
      "provenance": [],
      "toc_visible": true,
      "include_colab_link": true
    },
    "accelerator": "GPU"
  },
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "view-in-github",
        "colab_type": "text"
      },
      "source": [
        "<a href=\"https://colab.research.google.com/github/dvschultz/ml-art-colabs/blob/master/ML4A_TSNE.ipynb\" target=\"_parent\"><img src=\"https://colab.research.google.com/assets/colab-badge.svg\" alt=\"Open In Colab\"/></a>"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "vRrpCy8k6gKT"
      },
      "source": [
        "This notebook is modified from the [ML4A](https://github.com/ml4a/ml4a/blob/master/examples/info_retrieval/image-tsne.ipynb) repo to work with a custom dataset.\n",
        "\n",
        "Make sure you set `Runtime` > `Change runtime type` > `Hardware Accelerator` to GPU."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "ZlXNr2SB5HNq"
      },
      "source": [
        "# Image t-SNE\n",
        "\n",
        "This notebook will take you through the process of generating a [t-SNE](https://en.wikipedia.org/wiki/T-distributed_stochastic_neighbor_embedding) of a set of images, using a feature vector for each image derived from the activations of the last fully-connected layer in a pre-trained convolutional neural network (convnet)."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "x4WBZA140W2I"
      },
      "source": [
        "## Feature extraction and reverse image search\n",
        "This notebook will guide you through the procedure of analyzing a large set of images using a pre-trained convolutional network, extracting feature vectors for each one which represent each image.\n",
        "\n",
        "After the analysis is done, we will review some retrieval tasks that you can do with such an analysis. The main task will be that of \"reverse image search,\" which refers to searching for the most similar set of images to some query image."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "z_Nuzyg00s5F"
      },
      "source": [
        "### Prepare a dataset\n",
        "You will need a folder of images to analyze. There is no limit to the size of the dataset, but for good results, at least 1000 would be good.\n",
        "\n",
        "I have uploaded a dataset to my Google Drive and will access it by conncting Drive to Colab. Run the cell below and walk thru the steps to connect your Drive account.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "3n15n3131CD4"
      },
      "source": [
        "from google.colab import drive\n",
        "drive.mount('/content/drive')"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "lDdjVVto1pJq"
      },
      "source": [
        "Next I’ll copy my zip file out of Drive and into Colab. Then I’ll unzip it."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "PRpDIJrq1dWC"
      },
      "source": [
        "!cp /content/drive/MyDrive/recycling_center.zip /content/recycling_center.zip\n",
        "!unzip /content/recycling_center.zip"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "4QF-ZyKl18IY"
      },
      "source": [
        "Next we’ll load some libraries and packages needed"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "collapsed": true,
        "id": "o4X-VHIa5HNt"
      },
      "source": [
        "%matplotlib inline\n",
        "import os\n",
        "import random\n",
        "import numpy as np\n",
        "import json\n",
        "import matplotlib.pyplot as plt\n",
        "import pickle\n",
        "from matplotlib.pyplot import imshow\n",
        "from PIL import Image\n",
        "from sklearn.manifold import TSNE\n",
        "\n",
        "import keras\n",
        "from keras.preprocessing import image\n",
        "from keras.applications.imagenet_utils import decode_predictions, preprocess_input\n",
        "from keras.models import Model"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "kY80XsAr2Esw"
      },
      "source": [
        "OK now we’re mostly set up. We can start analyzing our features. Let’s load a pretrained model called VGG16 to analyze our images."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "p6wpet9r15vO"
      },
      "source": [
        "model = keras.applications.VGG16(weights='imagenet', include_top=True)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "NDEdKthi2SKf"
      },
      "source": [
        "Let's look at the model summary to see how it's structured."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "GfNHCFgR2TF5"
      },
      "source": [
        "model.summary()"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "G0ZpvIYy2ou0"
      },
      "source": [
        "### Forwarding an image through the network\n",
        "In order to input an image into the network, it has to be pre-processed into a feature vector of the correct size. To help us do this, we will create a function load_image(path) which will handle the usual pre-processing steps: load an image from our file system and turn it into an input vector of the correct dimensions, those expected by VGG16, namely a color image of size 224x224."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "j6aGeV8B2tkL"
      },
      "source": [
        "def load_image(path):\n",
        "    img = image.load_img(path, target_size=model.input_shape[1:3])\n",
        "    x = image.img_to_array(img)\n",
        "    x = np.expand_dims(x, axis=0)\n",
        "    x = preprocess_input(x)\n",
        "    return img, x"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "xOclLAeU3BCA"
      },
      "source": [
        "We'll load an image in from our dataset, and take a look at its data vector."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "hASTeM2W3AUf"
      },
      "source": [
        "img, x = load_image(\"/content/recycling_center/10715801114.jpg\")\n",
        "print(\"shape of x: \", x.shape)\n",
        "print(\"data type: \", x.dtype)\n",
        "plt.imshow(img)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "3YCBzFvX3QEU"
      },
      "source": [
        "The shape of the image is [1, 224, 224, 3]. The reason it has the extra first dimension with 1 element is that the network can take batches of images to process them all simultaneously. So for example, 10 images can be propagated through the network if x has a shape of [10, 224, 224, 3].\n",
        "\n",
        "Let's get class predictions from this model. We forward x through model and then use the built-in decode_predictions to look up the class names."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "NUlWgTkD3S84"
      },
      "source": [
        "# forward the image through the network\n",
        "predictions = model.predict(x)\n",
        "\n",
        "# print out the \n",
        "for _, pred, prob in decode_predictions(predictions)[0]:\n",
        "    print(\"predicted %s with probability %0.3f\" % (pred, prob))"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "dznx4q7z3Yk4"
      },
      "source": [
        "### Feature extraction\n",
        "What we have in the model variable is a highly effective image classifier trained on the ImageNet database. We expect that the classifier must form a very effective representation of the image in order to be able to classify it with such high accuracy. We can use this to our advantage by re-purposing this for another task.\n",
        "\n",
        "What we do is we copy the model, but remove the last layer (the classification layer), so that the final layer of the new network, called feat_extractor is the second 4096-neuron fully-connected layer, \"fc2 (Dense)\".\n",
        "\n",
        "The way we do this is by instantiating a new model called feature_extractor which takes a reference to the desired input and output layers in our VGG16 model. Thus, feature_extractor's output is the layer just before the classification, the last 4096-neuron fully connected layer. It looks like a copy, but internally, all Keras is doing is making a pointer to each of these layers and not actually copying anything. Thus, the output \"prediction\" from feat_extractor will just be the layer fc2 from model.\n",
        "\n",
        "If we run the summary() function again, we see that the architecture of feat_extractor is identical to the original model, except the last layer has been removed. We also know that not just the architecture is the same, but the two have the same weights as well."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "ddvUz-v13gg5"
      },
      "source": [
        "feat_extractor = Model(inputs=model.input, outputs=model.get_layer(\"fc2\").output)\n",
        "feat_extractor.summary()"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "7LPtor6C3liq"
      },
      "source": [
        "Now let's see the feature extractor in action. We pass the same image from before into it, and look at the results. The predict function returns an array with one element per image (in our case, there is just one). Each element contains a 4096-element array, which is the activations of the last fully-connected layer fc2 in VGG16. Let's plot the array as well."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "q2_lUpoH3qiQ"
      },
      "source": [
        "img, x = load_image(\"/content/recycling_center/101373297.jpg\")\n",
        "feat = feat_extractor.predict(x)\n",
        "\n",
        "plt.figure(figsize=(16,4))\n",
        "plt.plot(feat[0])"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "unTN7KI_3vy7"
      },
      "source": [
        "Our expectation is that the fc2 activations form a very good representation of the image, such that similar images should produce similar activations. In other words, the fc2 activations of two images which have similar content should be very close to each other. We can exploit this to do information retrieval.\n",
        "\n",
        "In the next cell, we will open a folder of images for analysis. First, the next cell will just recursively crawl the folder specified by image_path looking for images of the extensions inside of image_extensions and then limiting them to a random subset of maximum max_num_images images. Change these variables if you wish to change the target images."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "0FbaphwG-wuN"
      },
      "source": [
        "!rm /content/recycling_center/2111608980.jpg\n",
        "!rm /content/recycling_center/2086548971.jpg"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "U-txv3W63242"
      },
      "source": [
        "images_path = '/content/recycling_center'\n",
        "image_extensions = ['.jpg', '.png', '.jpeg']   # case-insensitive (upper/lower doesn't matter)\n",
        "max_num_images = 10000\n",
        "\n",
        "images = [os.path.join(dp, f) for dp, dn, filenames in os.walk(images_path) for f in filenames if os.path.splitext(f)[1].lower() in image_extensions]\n",
        "if max_num_images < len(images):\n",
        "    images = [images[i] for i in sorted(random.sample(xrange(len(images)), max_num_images))]\n",
        "\n",
        "print(\"keeping %d images to analyze\" % len(images))"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "raP5rrqu37kD"
      },
      "source": [
        "In the next cell, we will begin a loop which will open each image, extract its feature vector, and append it to a list called features which will contain our activations for each image. This process may take a long time depending on your graphics card, so you may need to leave it running for as much as a few hours. On a good graphics card, this process may only take a half hour or so. Every 1000 images, you will receive a notification print-out."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "fzifqIDK4CaW"
      },
      "source": [
        "import time\n",
        "tic = time.perf_counter()\n",
        "\n",
        "features = []\n",
        "for i, image_path in enumerate(images):\n",
        "    if i % 500 == 0:\n",
        "        toc = time.clock()\n",
        "        elap = toc-tic;\n",
        "        print(\"analyzing image %d / %d. Time: %4.4f seconds.\" % (i, len(images),elap))\n",
        "        tic = time.clock()\n",
        "    print(image_path)\n",
        "    img, x = load_image(image_path);\n",
        "    feat = feat_extractor.predict(x)[0]\n",
        "    features.append(feat)\n",
        "\n",
        "print('finished extracting features for %d images' % len(images))"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "6zPuNl4s4PC7"
      },
      "source": [
        "Alone, these activations provide a good representation, but it is a good idea to do one more step before using these as our feature vectors, which is to do a principal component analysis (PCA) to reduce the dimensionality of our feature vectors down to 300. We apply PCA for two reasons: 1) the 4096-bit feature vector may have some redundancy in it, such that multiple elements in the vector are highly correlated or similar. This would skew similarity comparisons towards those over-represented features. 2) Operating over 4096 elements is inefficient both in terms of space/memory requirements and processor speed, and it would be better for us if we can reduce the length of these vectors but maintain the same effective representation. PCA allows us to do this by reducing the dimensionality down of the feature vectors from 4096 to much less, but maintain a representation which is still faithful to the original data, by preserving the relative inter-point distance.\n",
        "\n",
        "Thus, PCA reduces the amount of redundancy in our features (from duplicate or highly-correlated features), speeds up computation over them, and reduces the amount of memory they take up.\n",
        "\n",
        "The next cell will instantiate a PCA object, which we will then fit our data to, choosing to keep the top 300 principal components. This may take a few minutes."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "IbEAbMRW4SA4"
      },
      "source": [
        "from sklearn.decomposition import PCA\n",
        "\n",
        "features = np.array(features)\n",
        "pca = PCA(n_components=300)\n",
        "pca.fit(features)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "tG6ScKUP4erx"
      },
      "source": [
        "The pca object stores the actual transformation matrix which was fit in the previous cell. We can now use it to transform any original feature vector (of length 4096) into a reduced 300-dimensional feature vector in the principal component space found by the PCA.\n",
        "\n",
        "So we take our original feature vectors, and transform them to the new space."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "69EVu0JU4fkM"
      },
      "source": [
        "pca_features = pca.transform(features)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "CD1yxi5X44Ha"
      },
      "source": [
        "We need to save both the image features matrix (the PCA-reduced features, not the originals), as well as the array containing the paths to each image, to make sure we can line up the images to their corresponding vectors. It is also best to save pca itself so we can project new images into the space if we wish. We can save everything to disk using pickle."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "w8mvWqS74sCe"
      },
      "source": [
        "import pickle\n",
        "\n",
        "pickle.dump([images, pca_features, pca], open('/content/features.pkl', 'wb'))"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "InzsY3BO5HNv"
      },
      "source": [
        "## TSNE\n",
        "\n",
        "First, we will load our image paths and feature vectors from the previous notebook step into memory. We can print their contents to get an idea of what they look like."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "Q8ajp-mt5HNw"
      },
      "source": [
        "images, pca_features, pca = pickle.load(open('/content/features.pkl', 'rb'))\n",
        "\n",
        "for img, f in list(zip(images, pca_features))[0:5]:\n",
        "    print(\"image: %s, features: %0.2f,%0.2f,%0.2f,%0.2f... \"%(img, f[0], f[1], f[2], f[3]))"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "5aaxEdu95HNx"
      },
      "source": [
        "Although in principle t-SNE works with any number of images, it's difficult to place that many tiles in a single image. So instead, we will take a random subset of 1000 images and plot those on a t-SNE instead. This step is optional, or you can try changing `num_images_to_plot`."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "QaNfNe7N5HNx"
      },
      "source": [
        "num_images_to_plot = 1000\n",
        "\n",
        "if len(images) > num_images_to_plot:\n",
        "    sort_order = sorted(random.sample(range(len(images)), num_images_to_plot))\n",
        "    images = [images[i] for i in sort_order]\n",
        "    pca_features = [pca_features[i] for i in sort_order]"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "CPhp6HC25HNy"
      },
      "source": [
        "It is usually a good idea to first run the vectors through a faster dimensionality reduction technique like [principal component analysis](https://en.wikipedia.org/wiki/Principal_component_analysis) to project your data into an intermediate lower-dimensional space before using t-SNE. This improves accuracy, and cuts down on runtime since PCA is more efficient than t-SNE. Since we have already projected our data down with PCA in the previous notebook, we can proceed straight to running the t-SNE on the feature vectors. Run the command in the following cell, taking note of the arguments:\n",
        "\n",
        "- `n_components` is the number of dimensions to project down to. In principle it can be anything, but in practice t-SNE is almost always used to project to 2 or 3 dimensions for visualization purposes.\n",
        "- `learning_rate` is the step size for iterations. You usually won't need to adjust this much, but your results may vary slightly. \n",
        "- `perplexity` refers to the number of independent clusters or zones t-SNE will attempt to fit points around. Again, it is relatively robust to large changes, and usually 20-50 works best. \n",
        "- `angle` controls the speed vs accuracy tradeoff. Lower angle means better accuracy but slower, although in practice, there is usually little improvement below a certain threshold."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "9GcSzKv25HNz"
      },
      "source": [
        "X = np.array(pca_features)\n",
        "tsne = TSNE(n_components=2, learning_rate=150, perplexity=30, angle=0.2, verbose=2).fit_transform(X)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "xwnnzjXw5HN0"
      },
      "source": [
        "Internally, t-SNE uses an iterative approach, making small (or sometimes large) adjustments to the points. By default, t-SNE will go a maximum of 1000 iterations, but in practice, it often terminates early because it has found a locally optimal (good enough) embedding.\n",
        "\n",
        "The variable `tsne` contains an array of unnormalized 2d points, corresponding to the embedding. In the next cell, we normalize the embedding so that lies entirely in the range (0,1)."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "collapsed": true,
        "id": "_-oGnNMH5HN1"
      },
      "source": [
        "tx, ty = tsne[:,0], tsne[:,1]\n",
        "tx = (tx-np.min(tx)) / (np.max(tx) - np.min(tx))\n",
        "ty = (ty-np.min(ty)) / (np.max(ty) - np.min(ty))"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "qtAji0Pe5HN2"
      },
      "source": [
        "Finally, we will compose a new RGB image where the set of images have been drawn according to the t-SNE results. Adjust `width` and `height` to set the size in pixels of the full image, and set `max_dim` to the pixel size (on the largest size) to scale images to."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "7DMwtgov5HN3"
      },
      "source": [
        "width = 4000\n",
        "height = 3000\n",
        "max_dim = 100\n",
        "\n",
        "full_image = Image.new('RGBA', (width, height))\n",
        "for img, x, y in zip(images, tx, ty):\n",
        "    tile = Image.open(img)\n",
        "    rs = max(1, tile.width/max_dim, tile.height/max_dim)\n",
        "    tile = tile.resize((int(tile.width/rs), int(tile.height/rs)), Image.ANTIALIAS)\n",
        "    full_image.paste(tile, (int((width-max_dim)*x), int((height-max_dim)*y)), mask=tile.convert('RGBA'))\n",
        "\n",
        "matplotlib.pyplot.figure(figsize = (16,12))\n",
        "imshow(full_image)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "TiSqGs9p5HN4"
      },
      "source": [
        "You can save the image to disk:"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "collapsed": true,
        "id": "n-z9x1rA5HN4"
      },
      "source": [
        "full_image.save(\"/content/recycling_center-tsne.png\")"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "nKtOzE_r5HN5"
      },
      "source": [
        "Now that we have generated our t-SNE, one more nice thing we can optionally do is to take the 2d embedding and assign it to a grid, using [RasterFairy](https://github.com/Quasimondo/RasterFairy). We can optionally choose a grid size of rows (`nx`) and columns (`ny`), which should be equal to the number of images you have. If it is less, then you can simply cut the `tsne` and `images` lists to be equal to `nx * ny`.\n",
        "\n",
        "If you omit the `target=(nx, ny)` argument, RasterFairy will automatically choose an optimal grid size to be as square-shaped as possible. RasterFairy also has options for embedding them in a grid with irregular borders as well (see the GitHub page for more details)."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "yYsSXrJU5HN5"
      },
      "source": [
        "You can also save the t-SNE points and their associated image paths for further processing in another environment."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "3bXQA_Xn5HN6"
      },
      "source": [
        "tsne_path = \"/content/example-tSNE-points.json\"\n",
        "\n",
        "data = [{\"path\":os.path.abspath(img), \"point\":[float(x), float(y)]} for img, x, y in zip(images, tx, ty)]\n",
        "with open(tsne_path, 'w') as outfile:\n",
        "    json.dump(data, outfile)\n",
        "\n",
        "print(\"saved t-SNE result to %s\" % tsne_path)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "vRLo1PkZE-_N"
      },
      "source": [
        "If you'd like to follow the last section which converts the t-SNE points to a grid assignment, you'll need [bmcfee's fork of Mario Klingemann's RasterFairy](https://github.com/bmcfee/RasterFairy), which can be installed with `pip` in with the following command."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "K_Ec5tN-E1pi"
      },
      "source": [
        "pip install -U git+https://github.com/bmcfee/RasterFairy/"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "lwVzzjx25HN6"
      },
      "source": [
        "import rasterfairy\n",
        "\n",
        "# nx * ny = 1000, the number of images\n",
        "nx = 40\n",
        "ny = 25\n",
        "\n",
        "# assign to grid\n",
        "grid_assignment = rasterfairy.transformPointCloud2D(tsne, target=(nx, ny))"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "GXtCnhty5HN7"
      },
      "source": [
        "Now finally, we can make a new image of our grid. Set the `tile_width` and `tile_height` variables according to how big you want the individual tile images to be. The resolution of the output image is `tile_width * nx` x `tile_height * ny`. The script will automatically center-crop all the tiles to match the aspect ratio of `tile_width / tile_height`."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "scrolled": false,
        "id": "4RGyLkyU5HN7"
      },
      "source": [
        "tile_width = 72\n",
        "tile_height = 56\n",
        "\n",
        "full_width = tile_width * nx\n",
        "full_height = tile_height * ny\n",
        "aspect_ratio = float(tile_width) / tile_height\n",
        "\n",
        "grid_image = Image.new('RGB', (full_width, full_height))\n",
        "\n",
        "for img, grid_pos in zip(images, grid_assignment[0]):\n",
        "    idx_x, idx_y = grid_pos\n",
        "    x, y = tile_width * idx_x, tile_height * idx_y\n",
        "    tile = Image.open(img)\n",
        "    tile_ar = float(tile.width) / tile.height  # center-crop the tile to match aspect_ratio\n",
        "    if (tile_ar > aspect_ratio):\n",
        "        margin = 0.5 * (tile.width - aspect_ratio * tile.height)\n",
        "        tile = tile.crop((margin, 0, margin + aspect_ratio * tile.height, tile.height))\n",
        "    else:\n",
        "        margin = 0.5 * (tile.height - float(tile.width) / aspect_ratio)\n",
        "        tile = tile.crop((0, margin, tile.width, margin + float(tile.width) / aspect_ratio))\n",
        "    tile = tile.resize((tile_width, tile_height), Image.ANTIALIAS)\n",
        "    grid_image.paste(tile, (int(x), int(y)))\n",
        "\n",
        "matplotlib.pyplot.figure(figsize = (16,12))\n",
        "imshow(grid_image)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "14qe0e1h5HN8"
      },
      "source": [
        "Finally, we can save the gridded t-SNE to disk as well."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "collapsed": true,
        "id": "Rl5lYMt-5HN9"
      },
      "source": [
        "grid_image.save(\"/content/recycling_center-tsne-2dgrid.jpg\")"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "Vlus7ltBFrj0"
      },
      "source": [
        ""
      ],
      "execution_count": null,
      "outputs": []
    }
  ]
}