tutorial_nbs/00_How_to_efficiently_work_with_very_large_numpy_arrays.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "created by Ignacio Oguiza - email: oguiza@timeseriesAI.co"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## How to efficiently work with (very large) Numpy Arrays? 👷‍♀️"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Sometimes we need to work with some very large numpy arrays that don't fit in memory. I'd like to share with you a way that works well for me."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Import libraries 📚"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# # **************** UNCOMMENT AND RUN THIS CELL IF YOU NEED TO INSTALL/ UPGRADE TSAI ****************\n",
    "# stable = True # Set to True for latest pip version or False for main branch in GitHub\n",
    "# !pip install {\"tsai -U\" if stable else \"git+https://github.com/timeseriesAI/tsai.git\"} >> /dev/null"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "os             : Darwin\n",
      "os version     : 19.6.0\n",
      "python         : 3.6.13\n",
      "tsai           : 0.2.20\n",
      "fastai         : 2.5.2\n",
      "fastcore       : 1.3.26\n",
      "numba          : 0.53.1\n",
      "torch          : 1.9.0\n",
      "n_cpus         : 8\n",
      "device         : cpu\n"
     ]
    }
   ],
   "source": [
    "from tsai.all import *\n",
    "my_setup()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Introduction 🤝"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "I normally work with time series data. I made the decision to use numpy arrays to store my data since they can easily handle multiple dimensions, and are really very efficient.\n",
    "\n",
    "But sometimes datasets are really big (many GBs) and don't fit in memory. So I started looking around and found something that works very well: [**np.memmap**](https://docs.scipy.org/doc/numpy/reference/generated/numpy.memmap.html). Conceptually they work as arrays on disk, and that's how I often call them.\n",
    "\n",
    "np.memmap creates a map to numpy arrays you have previously saved on disk, so that you can efficiently access small segments of those (small or large) files on disk, without reading the entire file into memory. And that's exactly what we need with deep learning, be able to quickly create a batch in memory, without reading the entire file (that is stored on disk). \n",
    "\n",
    "The best analogy I've found are image files. You may have a very large dataset on disk (that far exceeds your RAM). In order to create your DL datasets, what you pass are the paths to each individual file, so that you can then load a few images and create a batch on demand.\n",
    "\n",
    "You can view np.memmap as the path collection that can be used to load numpy data on demand when you need to create a batch.\n",
    "\n",
    "So let's see how you can work with larger than RAM arrays on disk."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "On my laptop I have only 8GB of RAM."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "I will try to demonstrate how you can handle a 10 GB numpy array dataset in an efficient way. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Create and save a larger-than-memory array 🥴"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "I will now create a large numpy array that doesn't fit in memory. \n",
    "Since I don't have enough RAM, I'll create an empty array on disk, and then load data in chunks that fit in memory.\n",
    "\n",
    "⚠️ If you want to to experiment with large datasets, you may uncomment and run this code. **It will create a ~10GB file on your disk** (we'll delete it at the end of this notebook).\n",
    "\n",
    "In my laptop it took me around **2 mins to create the data.**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "\n",
       "    <div>\n",
       "        <style>\n",
       "            /* Turns off some styling */\n",
       "            progress {\n",
       "                /* gets rid of default border in Firefox and Opera. */\n",
       "                border: none;\n",
       "                /* Needs to be in here for Safari polyfill so background images work as expected. */\n",
       "                background-size: auto;\n",
       "            }\n",
       "            .progress-bar-interrupted, .progress-bar-interrupted::-webkit-progress-bar {\n",
       "                background: #F44336;\n",
       "            }\n",
       "        </style>\n",
       "      <progress value='10' class='' max='10' style='width:300px; height:20px; vertical-align: middle;'></progress>\n",
       "      100.00% [10/10 02:16<00:00]\n",
       "    </div>\n",
       "    "
      ],
      "text/plain": [
       "<IPython.core.display.HTML object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# path = Path('./data')\n",
    "# X = create_empty_array((100_000, 50, 512), fname='X_on_disk', path=path, mode='r+')\n",
    "\n",
    "# chunksize = 10_000\n",
    "# pbar = progress_bar(range(math.ceil(len(X) / chunksize)))\n",
    "# start = 0\n",
    "# for i in pbar:\n",
    "#     end = start + chunksize\n",
    "#     X[start:end] = np.random.rand(chunksize, X.shape[-2], X.shape[-1])\n",
    "#     start = end\n",
    "\n",
    "# # I will create a smaller array. Sinc this fits in memory, I don't need to use a memmap\n",
    "# y_fn = path/'y_on_disk.npy'\n",
    "# y = np.random.randint(0, 10, X.shape[0])\n",
    "# labels = np.array(['a', 'b', 'c', 'd', 'e', 'f', 'g', 'h', 'i', 'j'])\n",
    "# np.save(y_fn, labels[y])\n",
    "\n",
    "# del X, y"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Ok. So let's check the size of these files if they were in memory."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "X array:  10240000128 bytes (9.54 GB)\n",
      "y array:       400128 bytes (0.00 GB)\n"
     ]
    }
   ],
   "source": [
    "print(f'X array: {os.path.getsize(\"./data/X_on_disk.npy\"):12} bytes ({bytes2str(os.path.getsize(\"./data/X_on_disk.npy\"))})')\n",
    "print(f'y array: {os.path.getsize(\"./data/y_on_disk.npy\"):12} bytes ({bytes2str(os.path.getsize(\"./data/y_on_disk.npy\"))})')\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Load an array on disk (np.memmap) 🧠"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Remember I only have an 8 GB RAM on this laptop, so I couldn't load these datasets in memory.\n",
    "\n",
    "☣️ Actually I accidentally loaded the \"X_on_disk.npy\" file, and my laptop crashed so I had to reboot it!\n",
    "\n",
    "So let's now load data as arrays on disk (np.memmap). The way to do it is super simple, and very efficient. You just do it as you would with a normal array, but add an mmap_mode.\n",
    "\n",
    "There are 4 modes: \n",
    "\n",
    "- ‘r’\tOpen existing file for reading only.\n",
    "- ‘r+’\tOpen existing file for reading and writing.\n",
    "- ‘w+’\tCreate or overwrite existing file for reading and writing.\n",
    "- ‘c’\tCopy-on-write: assignments affect data in memory, but changes are not saved to disk. The file on disk is read-only.\n",
    "\n",
    "I normally use mode 'c' since I want to be able to make changes to data in memory (transforms for example), without affecting data on disk (same approach as with image data). This is the same thing you do with image files on disk, that are just read, and then modified in memory, without changing the file on disk.\n",
    "\n",
    "But if you also want to be able to modify data on disk, you can load the array with mmap_mode='r+'."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "X_on_disk = np.load('./data/X_on_disk.npy', mmap_mode='c')\n",
    "y_on_disk = np.load('./data/y_on_disk.npy', mmap_mode='c')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Fast load**: it only takes a few ms to \"load\" a memory map to a 10 GB array on disk.\n",
    "\n",
    "In fact, the only thing that is loaded is a map to the array stored on disk. That's why it's so fast."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Arrays on disk: main features 📀"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Very limited RAM usage"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(100000, 50, 512) (100000,)\n"
     ]
    }
   ],
   "source": [
    "print(X_on_disk.shape, y_on_disk.shape)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "X array on disk:          152 bytes (0.000 GB)\n",
      "y array on disk:          120 bytes (0.000 GB)\n"
     ]
    }
   ],
   "source": [
    "print(f'X array on disk: {sys.getsizeof(X_on_disk):12} bytes ({bytes2str(sys.getsizeof(X_on_disk))})')\n",
    "print(f'y array on disk: {sys.getsizeof(y_on_disk):12} bytes ({bytes2str(sys.getsizeof(y_on_disk))})')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**152 bytes of RAM for a 10GB array**. This is the great benefit of arrays on disk.\n",
    "\n",
    "Arrays on disk barely use any RAM until they are accessed or sliced and an element is converted into a np.array or a tensor.\n",
    "\n",
    "This is equivalent to the size of file paths in images (very limited) compared to the files themselves (actual images). "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Types"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "np.memmap is a subclass of np.ndarray"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": null,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "isinstance(X_on_disk, np.ndarray)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "numpy.memmap"
      ]
     },
     "execution_count": null,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "type(X_on_disk)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Operations"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "With np.memmap you can perform the same operations you would with a normal numpy array. \n",
    "The most common operations you will perform in deep learning are:\n",
    "\n",
    "- slicing\n",
    "- calculating stats: mean and std\n",
    "- scaling (using normalize or standardize)\n",
    "- transformation into a tensor\n",
    "\n",
    "Once you get the slice from the array on disk, you'll convert it into a tensor, move to a GPU and performs operations there."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "⚠️ You need to be careful though not to convert the entire np.memmap to an array/ tensor if it's larger than your RAM. This will crash your computer unless you have enough RAM, so you would have to reboot!\n",
    "\n",
    "**DON'T DO THIS:  torch.from_numpy(X) or np.array(X)** unless you have enough RAM.\n",
    "\n",
    "To avoid issues during test, I created a smaller array on disk (that I can store in memory). When I want to test something I test it with that array first. It's important to always verify that the type output of your operations is np.memmap, which means data is still in memory."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Slicing"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To ensure you don't bring the entire array in memory (which may crash your computer) you can always work with slices of data, which is by the way how fastai works.\n",
    "\n",
    "If you use mode 'r' you can grab a sample and make changes to it, but this won't modify data on disk."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "memmap([[0.34803748, 0.04089686, 0.7991264 , ..., 0.37770525, 0.8188598 ,\n",
       "         0.17301634],\n",
       "        [0.6830333 , 0.6505463 , 0.37935442, ..., 0.8085311 , 0.04008082,\n",
       "         0.2903932 ],\n",
       "        [0.9397006 , 0.89364505, 0.16049282, ..., 0.15007935, 0.345074  ,\n",
       "         0.1874511 ],\n",
       "        ...,\n",
       "        [0.0391452 , 0.2724889 , 0.3282523 , ..., 0.9183814 , 0.67118037,\n",
       "         0.39444962],\n",
       "        [0.74068576, 0.16523811, 0.18562464, ..., 0.18377279, 0.3053631 ,\n",
       "         0.3580035 ],\n",
       "        [0.27491328, 0.7390123 , 0.9457232 , ..., 0.18275526, 0.6705215 ,\n",
       "         0.9529823 ]], dtype=float32)"
      ]
     },
     "execution_count": null,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "x = X_on_disk[0]\n",
    "x"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "It's important to note that **when we perform a math operation on a np.memmap (add, subtract, ...) the output is a np.array, and no longer a np.memmap.**\n",
    "\n",
    "⚠️ Remember you don't want to run this type of operations with a memmap larger than your RAM!! That's why I do it with a slice."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[1.3480375, 1.0408969, 1.7991264, ..., 1.3777052, 1.8188598,\n",
       "        1.1730163],\n",
       "       [1.6830332, 1.6505463, 1.3793545, ..., 1.808531 , 1.0400808,\n",
       "        1.2903932],\n",
       "       [1.9397006, 1.893645 , 1.1604928, ..., 1.1500794, 1.3450739,\n",
       "        1.1874511],\n",
       "       ...,\n",
       "       [1.0391452, 1.2724888, 1.3282523, ..., 1.9183815, 1.6711804,\n",
       "        1.3944496],\n",
       "       [1.7406857, 1.1652381, 1.1856246, ..., 1.1837728, 1.305363 ,\n",
       "        1.3580035],\n",
       "       [1.2749133, 1.7390122, 1.9457232, ..., 1.1827552, 1.6705215,\n",
       "        1.9529823]], dtype=float32)"
      ]
     },
     "execution_count": null,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "x = X_on_disk[0] + 1\n",
    "x"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "tensor([[1.3480, 1.0409, 1.7991,  ..., 1.3777, 1.8189, 1.1730],\n",
       "        [1.6830, 1.6505, 1.3794,  ..., 1.8085, 1.0401, 1.2904],\n",
       "        [1.9397, 1.8936, 1.1605,  ..., 1.1501, 1.3451, 1.1875],\n",
       "        ...,\n",
       "        [1.0391, 1.2725, 1.3283,  ..., 1.9184, 1.6712, 1.3944],\n",
       "        [1.7407, 1.1652, 1.1856,  ..., 1.1838, 1.3054, 1.3580],\n",
       "        [1.2749, 1.7390, 1.9457,  ..., 1.1828, 1.6705, 1.9530]])"
      ]
     },
     "execution_count": null,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "x = torch.from_numpy(X_on_disk[0])\n",
    "x2 = x + 1\n",
    "x2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As you can see, this doesn't affect the original np.memmap"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "memmap([[0.34803748, 0.04089686, 0.7991264 , ..., 0.37770525, 0.8188598 ,\n",
       "         0.17301634],\n",
       "        [0.6830333 , 0.6505463 , 0.37935442, ..., 0.8085311 , 0.04008082,\n",
       "         0.2903932 ],\n",
       "        [0.9397006 , 0.89364505, 0.16049282, ..., 0.15007935, 0.345074  ,\n",
       "         0.1874511 ],\n",
       "        ...,\n",
       "        [0.0391452 , 0.2724889 , 0.3282523 , ..., 0.9183814 , 0.67118037,\n",
       "         0.39444962],\n",
       "        [0.74068576, 0.16523811, 0.18562464, ..., 0.18377279, 0.3053631 ,\n",
       "         0.3580035 ],\n",
       "        [0.27491328, 0.7390123 , 0.9457232 , ..., 0.18275526, 0.6705215 ,\n",
       "         0.9529823 ]], dtype=float32)"
      ]
     },
     "execution_count": null,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "X_on_disk[0]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You can slice an array on disk by any axis, and it'll return a memmap. Slicing by any axis is very fast."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "memmap([[0.34803748, 0.04089686, 0.7991264 , ..., 0.37770525, 0.8188598 ,\n",
       "         0.17301634],\n",
       "        [0.6830333 , 0.6505463 , 0.37935442, ..., 0.8085311 , 0.04008082,\n",
       "         0.2903932 ],\n",
       "        [0.9397006 , 0.89364505, 0.16049282, ..., 0.15007935, 0.345074  ,\n",
       "         0.1874511 ],\n",
       "        ...,\n",
       "        [0.0391452 , 0.2724889 , 0.3282523 , ..., 0.9183814 , 0.67118037,\n",
       "         0.39444962],\n",
       "        [0.74068576, 0.16523811, 0.18562464, ..., 0.18377279, 0.3053631 ,\n",
       "         0.3580035 ],\n",
       "        [0.27491328, 0.7390123 , 0.9457232 , ..., 0.18275526, 0.6705215 ,\n",
       "         0.9529823 ]], dtype=float32)"
      ]
     },
     "execution_count": null,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "X_on_disk[0]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "memmap([[0.34803748, 0.04089686, 0.7991264 , ..., 0.37770525, 0.8188598 ,\n",
       "         0.17301634],\n",
       "        [0.9259249 , 0.514623  , 0.50216776, ..., 0.8823718 , 0.561646  ,\n",
       "         0.25591376],\n",
       "        [0.06298492, 0.10742943, 0.43376994, ..., 0.01061168, 0.3993792 ,\n",
       "         0.5877482 ],\n",
       "        ...,\n",
       "        [0.36538476, 0.6251516 , 0.13214637, ..., 0.56368643, 0.03602772,\n",
       "         0.02040654],\n",
       "        [0.7697917 , 0.06593986, 0.12318378, ..., 0.24622898, 0.4352764 ,\n",
       "         0.8795757 ],\n",
       "        [0.30351886, 0.05458342, 0.18446152, ..., 0.00465104, 0.35671628,\n",
       "         0.12464925]], dtype=float32)"
      ]
     },
     "execution_count": null,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "X_on_disk[:, 0]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "However, bear in mind that if you use multiple indices, the output will be a regular numpy array. This is important as it will use more RAM. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[[0.34803748, 0.04089686, 0.7991264 , ..., 0.37770525,\n",
       "         0.8188598 , 0.17301634],\n",
       "        [0.6830333 , 0.6505463 , 0.37935442, ..., 0.8085311 ,\n",
       "         0.04008082, 0.2903932 ],\n",
       "        [0.9397006 , 0.89364505, 0.16049282, ..., 0.15007935,\n",
       "         0.345074  , 0.1874511 ],\n",
       "        ...,\n",
       "        [0.0391452 , 0.2724889 , 0.3282523 , ..., 0.9183814 ,\n",
       "         0.67118037, 0.39444962],\n",
       "        [0.74068576, 0.16523811, 0.18562464, ..., 0.18377279,\n",
       "         0.3053631 , 0.3580035 ],\n",
       "        [0.27491328, 0.7390123 , 0.9457232 , ..., 0.18275526,\n",
       "         0.6705215 , 0.9529823 ]],\n",
       "\n",
       "       [[0.9259249 , 0.514623  , 0.50216776, ..., 0.8823718 ,\n",
       "         0.561646  , 0.25591376],\n",
       "        [0.81914663, 0.79375005, 0.65909016, ..., 0.884909  ,\n",
       "         0.23646063, 0.5160194 ],\n",
       "        [0.99880844, 0.6775859 , 0.16700691, ..., 0.84936655,\n",
       "         0.051814  , 0.20492136],\n",
       "        ...,\n",
       "        [0.4639562 , 0.41425797, 0.49373862, ..., 0.58005303,\n",
       "         0.17869665, 0.97369766],\n",
       "        [0.56871843, 0.20894419, 0.6577927 , ..., 0.8131621 ,\n",
       "         0.88691944, 0.5066556 ],\n",
       "        [0.3088769 , 0.883312  , 0.8087257 , ..., 0.9424154 ,\n",
       "         0.26620054, 0.32207966]]], dtype=float32)"
      ]
     },
     "execution_count": null,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "X_on_disk[[0,1]]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Unless you use a slice with consecutive indices like this:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "memmap([[[0.34803748, 0.04089686, 0.7991264 , ..., 0.37770525,\n",
       "          0.8188598 , 0.17301634],\n",
       "         [0.6830333 , 0.6505463 , 0.37935442, ..., 0.8085311 ,\n",
       "          0.04008082, 0.2903932 ],\n",
       "         [0.9397006 , 0.89364505, 0.16049282, ..., 0.15007935,\n",
       "          0.345074  , 0.1874511 ],\n",
       "         ...,\n",
       "         [0.0391452 , 0.2724889 , 0.3282523 , ..., 0.9183814 ,\n",
       "          0.67118037, 0.39444962],\n",
       "         [0.74068576, 0.16523811, 0.18562464, ..., 0.18377279,\n",
       "          0.3053631 , 0.3580035 ],\n",
       "         [0.27491328, 0.7390123 , 0.9457232 , ..., 0.18275526,\n",
       "          0.6705215 , 0.9529823 ]],\n",
       "\n",
       "        [[0.9259249 , 0.514623  , 0.50216776, ..., 0.8823718 ,\n",
       "          0.561646  , 0.25591376],\n",
       "         [0.81914663, 0.79375005, 0.65909016, ..., 0.884909  ,\n",
       "          0.23646063, 0.5160194 ],\n",
       "         [0.99880844, 0.6775859 , 0.16700691, ..., 0.84936655,\n",
       "          0.051814  , 0.20492136],\n",
       "         ...,\n",
       "         [0.4639562 , 0.41425797, 0.49373862, ..., 0.58005303,\n",
       "          0.17869665, 0.97369766],\n",
       "         [0.56871843, 0.20894419, 0.6577927 , ..., 0.8131621 ,\n",
       "          0.88691944, 0.5066556 ],\n",
       "         [0.3088769 , 0.883312  , 0.8087257 , ..., 0.9424154 ,\n",
       "          0.26620054, 0.32207966]]], dtype=float32)"
      ]
     },
     "execution_count": null,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "X_on_disk[:2]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This continues to be a memmap"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "There's a trick we can use avoid this making use of the excellent new L class in fastai. It is to **itemify** the np.memmap/s. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def itemify(*x): return L(*x).zip()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To itemify one or several np.memmap/s is very fast. Let's see how long it takes with a 10 GB array."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "X_on_disk_as_items = itemify(X_on_disk)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "5 seconds to return individual records on disk! Bear in mind you only need to perform this once!\n",
    "\n",
    "So now, you can select multiple items at the same time, and they will all still be on disk:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(#2) [(memmap([[0.34803748, 0.04089686, 0.7991264 , ..., 0.37770525, 0.8188598 ,\n",
       "         0.17301634],\n",
       "        [0.6830333 , 0.6505463 , 0.37935442, ..., 0.8085311 , 0.04008082,\n",
       "         0.2903932 ],\n",
       "        [0.9397006 , 0.89364505, 0.16049282, ..., 0.15007935, 0.345074  ,\n",
       "         0.1874511 ],\n",
       "        ...,\n",
       "        [0.0391452 , 0.2724889 , 0.3282523 , ..., 0.9183814 , 0.67118037,\n",
       "         0.39444962],\n",
       "        [0.74068576, 0.16523811, 0.18562464, ..., 0.18377279, 0.3053631 ,\n",
       "         0.3580035 ],\n",
       "        [0.27491328, 0.7390123 , 0.9457232 , ..., 0.18275526, 0.6705215 ,\n",
       "         0.9529823 ]], dtype=float32),),(memmap([[0.9259249 , 0.514623  , 0.50216776, ..., 0.8823718 , 0.561646  ,\n",
       "         0.25591376],\n",
       "        [0.81914663, 0.79375005, 0.65909016, ..., 0.884909  , 0.23646063,\n",
       "         0.5160194 ],\n",
       "        [0.99880844, 0.6775859 , 0.16700691, ..., 0.84936655, 0.051814  ,\n",
       "         0.20492136],\n",
       "        ...,\n",
       "        [0.4639562 , 0.41425797, 0.49373862, ..., 0.58005303, 0.17869665,\n",
       "         0.97369766],\n",
       "        [0.56871843, 0.20894419, 0.6577927 , ..., 0.8131621 , 0.88691944,\n",
       "         0.5066556 ],\n",
       "        [0.3088769 , 0.883312  , 0.8087257 , ..., 0.9424154 , 0.26620054,\n",
       "         0.32207966]], dtype=float32),)]"
      ]
     },
     "execution_count": null,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "X_on_disk_as_items[0,1]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You can also itemify several items at once: X and y for example. When you slice the list, you'll get tuples."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "Xy_on_disk_as_items = itemify(X_on_disk, y_on_disk)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(#2) [(memmap([[0.34803748, 0.04089686, 0.7991264 , ..., 0.37770525, 0.8188598 ,\n",
       "         0.17301634],\n",
       "        [0.6830333 , 0.6505463 , 0.37935442, ..., 0.8085311 , 0.04008082,\n",
       "         0.2903932 ],\n",
       "        [0.9397006 , 0.89364505, 0.16049282, ..., 0.15007935, 0.345074  ,\n",
       "         0.1874511 ],\n",
       "        ...,\n",
       "        [0.0391452 , 0.2724889 , 0.3282523 , ..., 0.9183814 , 0.67118037,\n",
       "         0.39444962],\n",
       "        [0.74068576, 0.16523811, 0.18562464, ..., 0.18377279, 0.3053631 ,\n",
       "         0.3580035 ],\n",
       "        [0.27491328, 0.7390123 , 0.9457232 , ..., 0.18275526, 0.6705215 ,\n",
       "         0.9529823 ]], dtype=float32), 'h'),(memmap([[0.9259249 , 0.514623  , 0.50216776, ..., 0.8823718 , 0.561646  ,\n",
       "         0.25591376],\n",
       "        [0.81914663, 0.79375005, 0.65909016, ..., 0.884909  , 0.23646063,\n",
       "         0.5160194 ],\n",
       "        [0.99880844, 0.6775859 , 0.16700691, ..., 0.84936655, 0.051814  ,\n",
       "         0.20492136],\n",
       "        ...,\n",
       "        [0.4639562 , 0.41425797, 0.49373862, ..., 0.58005303, 0.17869665,\n",
       "         0.97369766],\n",
       "        [0.56871843, 0.20894419, 0.6577927 , ..., 0.8131621 , 0.88691944,\n",
       "         0.5066556 ],\n",
       "        [0.3088769 , 0.883312  , 0.8087257 , ..., 0.9424154 , 0.26620054,\n",
       "         0.32207966]], dtype=float32), 'i')]"
      ]
     },
     "execution_count": null,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Xy_on_disk_as_items[0, 1]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Slicing is very fast, even if there are 100.000 samples."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "The slowest run took 21.11 times longer than the fastest. This could mean that an intermediate result is being cached.\n",
      "1000000 loops, best of 3: 1.93 µs per loop\n"
     ]
    }
   ],
   "source": [
    "# axis 0\n",
    "%timeit X_on_disk[0]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "The slowest run took 16.79 times longer than the fastest. This could mean that an intermediate result is being cached.\n",
      "100000 loops, best of 3: 1.96 µs per loop\n"
     ]
    }
   ],
   "source": [
    "# axis 1\n",
    "%timeit X_on_disk[:, 0]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "The slowest run took 19.08 times longer than the fastest. This could mean that an intermediate result is being cached.\n",
      "1000000 loops, best of 3: 1.94 µs per loop\n"
     ]
    }
   ],
   "source": [
    "# axis 2\n",
    "%timeit X_on_disk[..., 0]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "The slowest run took 20.50 times longer than the fastest. This could mean that an intermediate result is being cached.\n",
      "1000000 loops, best of 3: 1.97 µs per loop\n"
     ]
    }
   ],
   "source": [
    "# aixs 0,1\n",
    "%timeit X_on_disk[0, 0]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To compare how fast you can slice a np.memmap, let's create a smaller array that I can fit in memory (X_in_memory). This is 10 times smaller (100 MB) than the one on disk."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "X_in_memory_small = np.random.rand(10000, 50, 512)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "The slowest run took 46.56 times longer than the fastest. This could mean that an intermediate result is being cached.\n",
      "10000000 loops, best of 3: 183 ns per loop\n"
     ]
    }
   ],
   "source": [
    "%timeit X_in_memory_small[0]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's create the same array on disk. It's super simple:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "np.save('./data/X_on_disk_small.npy', X_in_memory_small)\n",
    "X_on_disk_small = np.load('./data/X_on_disk_small.npy', mmap_mode='c')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "The slowest run took 20.41 times longer than the fastest. This could mean that an intermediate result is being cached.\n",
      "1000000 loops, best of 3: 1.95 µs per loop\n"
     ]
    }
   ],
   "source": [
    "%timeit X_on_disk_small[0]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This is approximately 10 times slower than having arrays on disk, although it's still pretty fast.\n",
    "\n",
    "However, if we use the itemified version, it's much faster:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "The slowest run took 15.68 times longer than the fastest. This could mean that an intermediate result is being cached.\n",
      "1000000 loops, best of 3: 877 ns per loop\n"
     ]
    }
   ],
   "source": [
    "%timeit X_on_disk_as_items[0]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This is much better! So now you can access 1 of multiple items on disk with a pretty good performance."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Calculating stats: mean and std"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Another benefit of using arrays on disk is that you can calculate the mean and std deviation of the entire dataset. \n",
    "\n",
    "It takes a considerable time since the array is very big (10GB), but it's feasible:\n",
    "\n",
    "- mean (0.4999966):  1 min 45 s\n",
    "- std  (0.2886839): 11 min 43 s \n",
    "\n",
    "in my laptop. \n",
    "If you need them, you could calculate these stats once, and store the results (similar to ImageNet stats).\n",
    "However, you usually need to claculate these metrics for labeled (train) datasets, that tend to be smaller."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# X_mean = X_on_disk.mean()\n",
    "# X_mean"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# X_std = X_on_disk.std()\n",
    "# X_std"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Conversion into a tensor"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Conversion from an array on disk slice into a tensor is also very fast:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "tensor([[0.3480, 0.0409, 0.7991,  ..., 0.3777, 0.8189, 0.1730],\n",
       "        [0.6830, 0.6505, 0.3794,  ..., 0.8085, 0.0401, 0.2904],\n",
       "        [0.9397, 0.8936, 0.1605,  ..., 0.1501, 0.3451, 0.1875],\n",
       "        ...,\n",
       "        [0.0391, 0.2725, 0.3283,  ..., 0.9184, 0.6712, 0.3944],\n",
       "        [0.7407, 0.1652, 0.1856,  ..., 0.1838, 0.3054, 0.3580],\n",
       "        [0.2749, 0.7390, 0.9457,  ..., 0.1828, 0.6705, 0.9530]])"
      ]
     },
     "execution_count": null,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "torch.from_numpy(X_on_disk[0])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "X_on_disk_small_0 = X_on_disk_small[0]\n",
    "X_in_memory_small_0 = X_in_memory_small[0]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "The slowest run took 44.11 times longer than the fastest. This could mean that an intermediate result is being cached.\n",
      "1000000 loops, best of 3: 1.13 µs per loop\n"
     ]
    }
   ],
   "source": [
    "%timeit torch.from_numpy(X_on_disk_small_0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "The slowest run took 46.27 times longer than the fastest. This could mean that an intermediate result is being cached.\n",
      "1000000 loops, best of 3: 1.16 µs per loop\n"
     ]
    }
   ],
   "source": [
    "%timeit torch.from_numpy(X_in_memory_small_0 )"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "So it takes the same time to convert from numpy.memmap or from a np.array in memory."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Combined operations: slicing plus conversion to tensor"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's now check performance of the combined process: slicing plus conversion to a tensor. Based on what we've seen there are 3 options: \n",
    "\n",
    "- slice np.array in memory + conversion to tensor\n",
    "- slice np.memmap on disk + conversion to tensor\n",
    "- slice itemified np.memmap + converion to tensor"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "The slowest run took 76.50 times longer than the fastest. This could mean that an intermediate result is being cached.\n",
      "1000000 loops, best of 3: 1.44 µs per loop\n"
     ]
    }
   ],
   "source": [
    "%timeit torch.from_numpy(X_in_memory_small[0])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "The slowest run took 24.08 times longer than the fastest. This could mean that an intermediate result is being cached.\n",
      "100000 loops, best of 3: 3.89 µs per loop\n"
     ]
    }
   ],
   "source": [
    "%timeit torch.from_numpy(X_on_disk_small[0])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "X_on_disk_small_as_items = itemify(X_on_disk_small)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "The slowest run took 23.52 times longer than the fastest. This could mean that an intermediate result is being cached.\n",
      "100000 loops, best of 3: 2.35 µs per loop\n"
     ]
    }
   ],
   "source": [
    "%timeit torch.from_numpy(X_on_disk_small_as_items[0][0])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "So this last method is **almost as fast as having the array in memory**!! This is an excellent outcome, since slicing arrays in memory is a highly optimized operation. \n",
    "\n",
    "And we have the benefit of having access to very large datasets if needed."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Remove the arrays on disk"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Don't forget to remove the arrays you have created on disk."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "os.remove('./data/X_on_disk.npy')\n",
    "os.remove('./data/X_on_disk_small.npy')\n",
    "os.remove('./data/y_on_disk.npy')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Summary ✅"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We now have a very efficient way to work with very large numpy arrays.\n",
    "\n",
    "The process is very simple:\n",
    "\n",
    "- create and save the array on disk (as described before)\n",
    "- load it with a mmap_mode='c' if you want to be able to modify data in memory but not on dis, or 'r+ if you want to modify data both in memory and on disk.\n",
    "\n",
    "So my recommendation would be:\n",
    "\n",
    "- use numpy arrays in memory when possible (if your data fits in memory)\n",
    "- use numpy memmap (arrays on disk) when data doesn't fit. You will still have a great performance."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "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.9.12"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}