02.02 Universal functions, aggregation and broadcasting in NumPy.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Векторизованные вычисления\n",
    "До сих пор мы рассмотрели самые простые возможности NumPy. Теперь мы перейдем к наиболее важной особенности массивов NumPy - векторизованным вычислениям.\n",
    "\n",
    "Обычно для обработки элементов массива используются циклы. Рассмотрим пример вычисления обратного значения ($1/x$) для каждого элемента массива"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2017-12-05T04:39:22.517079Z",
     "start_time": "2017-12-05T04:39:22.117056Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([ 0.14285714,  0.25      ,  0.125     ,  0.2       ,  0.14285714])"
      ]
     },
     "execution_count": 1,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "import numpy as np\n",
    "np.random.seed(42)\n",
    "\n",
    "def reciprocals(values):\n",
    "    x = np.empty(len(values))\n",
    "    for i in range(len(values)):\n",
    "        x[i] = 1.0 / values[i]\n",
    "    return x\n",
    "\n",
    "values = np.random.randint(1, 10, size=5)\n",
    "reciprocals(values)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Давайте измерим время работы данной функции для очень большого массива с помощью magic команды %magic"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2017-12-05T04:39:59.355186Z",
     "start_time": "2017-12-05T04:39:22.521079Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "4.61 s ± 351 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)\n"
     ]
    }
   ],
   "source": [
    "big_array = np.random.randint(1, 1000, size=1000000)\n",
    "%timeit reciprocals(big_array)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Наименьшее время на вызов функции `reciprocals` составило ~2.5 секунды (у вас это время может отличаться). Это крайне много по сравнению, например, с языком C. Причина такого медленного выполнения в том, что при каждой итерации Python проверяет тип элемента массива и выбирает функцию, которая соответствует для данного типа оператора. Если тип каждого элемента был бы известен заранее, то таких накладных расходов не было.\n",
    "\n",
    "Вместо использования такого явного цикла NumPy позволяет использовать многие операторы напрямую с массивом. При этом оператор применится для каждого элемента массива. Такие вычисления называются векторизованными, так как мы применяем их сразу к целому вектору (т.е. массиву)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2017-12-05T04:39:59.367187Z",
     "start_time": "2017-12-05T04:39:59.358186Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[ 0.14285714  0.25        0.125       0.2         0.14285714]\n",
      "[ 0.14285714  0.25        0.125       0.2         0.14285714]\n"
     ]
    }
   ],
   "source": [
    "print(reciprocals(values))\n",
    "print(1.0 / values)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Как видно мы получили одинаковый результат для обоих вызовов. Давайте проверим время выполнения "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2017-12-05T04:40:04.344472Z",
     "start_time": "2017-12-05T04:39:59.372187Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "6.03 ms ± 605 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)\n"
     ]
    }
   ],
   "source": [
    "%timeit (1.0 / big_array)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Наименьшее время в ~830 раз меньше, чем при вызове функции `reciprocals` (~2.5 секунды против ~3 миллисекунд. Ваш результат может отличаться, но разница будет в таких же масштабах). Как видно, разница колоссальная.\n",
    "\n",
    "Векторизованные вычисления можно выполнять между двумя массивами. В примере ниже числитель и знаменатель дроби являются массивами"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2017-12-05T04:40:04.363473Z",
     "start_time": "2017-12-05T04:40:04.347472Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([ 6.        ,  2.5       ,  1.33333333,  0.75      ,  0.4       ])"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.arange(6, 1, -1) / np.arange(1, 6)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "В векторизованных вичислениях также могут участвовать многомерные массивы"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2017-12-05T04:40:04.392474Z",
     "start_time": "2017-12-05T04:40:04.372473Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[  1,   2,   4],\n",
       "       [  8,  16,  32],\n",
       "       [ 64, 128, 256]], dtype=int32)"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "x = np.arange(9).reshape((3, 3))\n",
    "2 ** x"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Можно строить выражения любой сложности для векторизованных вычислений. Давайте рассмотрим несколько примеров"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2017-12-05T04:40:04.415476Z",
     "start_time": "2017-12-05T04:40:04.398475Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([1, 2, 3, 4])"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "x = np.arange(1, 5)\n",
    "x"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2017-12-05T04:40:04.434477Z",
     "start_time": "2017-12-05T04:40:04.424476Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([4, 5, 6, 7])"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "x + 3"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2017-12-05T04:40:04.453478Z",
     "start_time": "2017-12-05T04:40:04.439477Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([3, 5, 7, 9])"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "x * 2 + 1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2017-12-05T04:40:04.474479Z",
     "start_time": "2017-12-05T04:40:04.457478Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([ 0.04761905,  0.0297619 ,  0.02380952,  0.02083333])"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "((x + 3) / 42) / (2 * x)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "NumPy включает в себя большое количество полезных функций, которые поддерживают векторизованные операции. Например, тригонометрические функции"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2017-12-05T04:40:04.491480Z",
     "start_time": "2017-12-05T04:40:04.481479Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "cos^2(x) + sin^2(x) = [ 1.  1.  1.  1.]\n"
     ]
    }
   ],
   "source": [
    "print('cos^2(x) + sin^2(x) =', np.cos(x) ** 2 + np.sin(x) ** 2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Или экспонента и логарифмы"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2017-12-05T04:40:04.520482Z",
     "start_time": "2017-12-05T04:40:04.499480Z"
    },
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "e^x        = [  2.71828183   7.3890561   20.08553692  54.59815003]\n",
      "e^x - 1    = [  1.71828183   6.3890561   19.08553692  53.59815003]\n",
      "ln(x)      = [ 0.          0.69314718  1.09861229  1.38629436]\n",
      "ln(1 + x)  = [ 0.69314718  1.09861229  1.38629436  1.60943791]\n",
      "log2(x)    = [ 0.         1.         1.5849625  2.       ]\n",
      "log10(x)   = [ 0.         1.         1.5849625  2.       ]\n"
     ]
    }
   ],
   "source": [
    "print('e^x        =', np.exp(x))\n",
    "print('e^x - 1    =', np.expm1(x))\n",
    "print('ln(x)      =', np.log(x))\n",
    "print('ln(1 + x)  =', np.log1p(x))\n",
    "print('log2(x)    =', np.log2(x))\n",
    "print('log10(x)   =', np.log2(x))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Список всех поддерживаемых функций можно посмотреть в документации NumPy. Некоторые специальные математические функции (для оптимизации, диффиренцирования, интегрирования и т.д.) также можно найти в модуле [SciPy](https://www.scipy.org/scipylib/index.html). Эти функции также поддерживают векторизованные вычисления.\n",
    "\n",
    "## Агрегирующие функции\n",
    "NumPy поддерживает простые агрегирующие функции такие как минимум, максимум и сумма"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2017-12-05T04:40:04.599486Z",
     "start_time": "2017-12-05T04:40:04.523482Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(1000000,)"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "x = np.random.normal(0, 100, size=1000000)\n",
    "x.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2017-12-05T04:40:04.622487Z",
     "start_time": "2017-12-05T04:40:04.602486Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "min: -480.871903296\n",
      "max: 462.132900248\n",
      "sum: -37376.6084239\n"
     ]
    }
   ],
   "source": [
    "print('min:', np.min(x))\n",
    "print('max:', np.max(x))\n",
    "print('sum:', np.sum(x))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Часто необходимо получить индекс минимального или максимального элемента"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2017-12-05T04:40:04.640488Z",
     "start_time": "2017-12-05T04:40:04.625488Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "argmin: 159908 value: -480.871903296\n",
      "argmax: 193533 value: 462.132900248\n"
     ]
    }
   ],
   "source": [
    "print('argmin:', np.argmin(x), 'value:', x[np.argmin(x)])\n",
    "print('argmax:', np.argmax(x), 'value:', x[np.argmax(x)])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Также есть набор стандартных статистик"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2017-12-05T04:40:04.699492Z",
     "start_time": "2017-12-05T04:40:04.643489Z"
    },
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "mean  :  -0.0373766084239\n",
      "median:  0.0814015827587\n",
      "std   :  100.004072348\n",
      "var   :  10000.8144862\n"
     ]
    }
   ],
   "source": [
    "print('mean  : ', np.mean(x))\n",
    "print('median: ', np.median(x))\n",
    "print('std   : ', np.std(x))\n",
    "print('var   : ', np.var(x))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Агрегирование многомерных массивов\n",
    "При агрегировании многомерных массивов можно указать ось (`axis`), по которому необходимо выполнить агрегацию"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2017-12-05T04:40:04.714493Z",
     "start_time": "2017-12-05T04:40:04.704492Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[6 4 5]\n",
      " [4 8 7]]\n"
     ]
    }
   ],
   "source": [
    "X = np.random.randint(0, 10, size=(2, 3))\n",
    "print(X)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2017-12-05T04:40:04.744494Z",
     "start_time": "2017-12-05T04:40:04.720493Z"
    },
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "total sum : 34\n",
      "column sum: [10 12 12]\n",
      "row sum   : [15 19]\n"
     ]
    }
   ],
   "source": [
    "print('total sum :', np.sum(X))\n",
    "print('column sum:', np.sum(X, axis=0))\n",
    "print('row sum   :', np.sum(X, axis=1))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2017-12-05T04:40:04.777496Z",
     "start_time": "2017-12-05T04:40:04.755495Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "maximum   : 8\n",
      "column max: [6 8 7]\n",
      "row max   : [6 8]\n"
     ]
    }
   ],
   "source": [
    "print('maximum   :', np.max(X))\n",
    "print('column max:', np.max(X, axis=0))\n",
    "print('row max   :', np.max(X, axis=1))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Список агригирующих функций\n",
    "\n",
    "|Функция            |   Без учета NaN     | Описание                                      |\n",
    "|-------------------|---------------------|-----------------------------------------------|\n",
    "| `np.min`          | `np.nanmin`         | Минимум                                       |\n",
    "| `np.max`          | `np.nanmax`         | Максимум                                      |\n",
    "| `np.argmin`       | `np.nanargmin`      | Индекс миниимального элемента массива         |\n",
    "| `np.argmax`       | `np.nanargmax`      | Индекс максимального элемента массива         |\n",
    "| `np.sum`          | `np.nansum`         | Сумма                                         |\n",
    "| `np.prod`         | `np.nanprod`        | Произведение                                  |\n",
    "| `np.mean`         | `np.nanmean`        | Среднее значение                              |\n",
    "| `np.median`       | `np.nanmedian`      | Медиана                                       |\n",
    "| `np.std`          | `np.nanstd`         | Стандартное отклонение                        |\n",
    "| `np.var`          | `np.nanvar`         | Дисперсия                                     |\n",
    "| `np.percentile`   | `np.nanpercentile`  | Перцентиль                                    |\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Кумулятивные функции\n",
    "В NumPy также есть кумулятивные функции, которые вычисляют накапливающее значение функции. Например, `np.cumsum` и `np.cumprod` вычисляют накапливающиеся сумму и произведение соответственно."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2017-12-05T04:40:04.793497Z",
     "start_time": "2017-12-05T04:40:04.781497Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "cumsum : [ 1  3  6 10]\n",
      "cumprod: [ 1  2  6 24]\n"
     ]
    }
   ],
   "source": [
    "x = np.arange(1, 5)\n",
    "print('cumsum :', np.cumsum(x))\n",
    "print('cumprod:', np.cumprod(x))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Ниже приведено вычисление для выражения $\\sum_{k=0}^{10} \\frac{1}{k!}$. Знаете что вычисляет данное выражение?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2017-12-05T04:40:04.810498Z",
     "start_time": "2017-12-05T04:40:04.796497Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "2.7182815255731922"
      ]
     },
     "execution_count": 21,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "x = np.arange(1, 10)\n",
    "1 + np.sum(1 / np.cumprod(x))"
   ]
  }
 ],
 "metadata": {
  "anaconda-cloud": {},
  "hide_input": false,
  "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.6.1"
  },
  "varInspector": {
   "cols": {
    "lenName": 16,
    "lenType": 16,
    "lenVar": 40
   },
   "kernels_config": {
    "python": {
     "delete_cmd_postfix": "",
     "delete_cmd_prefix": "del ",
     "library": "var_list.py",
     "varRefreshCmd": "print(var_dic_list())"
    },
    "r": {
     "delete_cmd_postfix": ") ",
     "delete_cmd_prefix": "rm(",
     "library": "var_list.r",
     "varRefreshCmd": "cat(var_dic_list()) "
    }
   },
   "types_to_exclude": [
    "module",
    "function",
    "builtin_function_or_method",
    "instance",
    "_Feature"
   ],
   "window_display": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}