examples/user_guide/15-Large_Data.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Working with large data using Datashader\n",
    "\n",
    "The various plotting-library backends supported by HoloViews, such as Matplotlib and Bokeh, each have a variety of limitations on the amount of data that is practical to work with.  Bokeh in particular mirrors your data directly into an HTML page viewable in your browser, which can cause problems when data sizes approach the limited memory available for each web page in current browsers.\n",
    "\n",
    "Luckily, a visualization of even the largest dataset will be constrained by the resolution of your display device, and so one approach to handling such data is to pre-render or rasterize the data into a fixed-size array or image *before* sending it to the backend.  The [Datashader](https://github.com/bokeh/datashader) library provides a high-performance big-data server-side rasterization pipeline that works seamlessly with HoloViews to support datasets that are orders of magnitude larger than those supported natively by the plotting-library backends, including millions or billions of points even on ordinary laptops.\n",
    "\n",
    "Here, we will see how and when to use Datashader with HoloViews Elements and Containers. For simplicity in this discussion we'll focus on simple synthetic datasets, but the [Datashader docs](http://datashader.org/topics) include a wide variety of real datasets that give a much better idea of the power of using Datashader with HoloViews, and [PyViz.org](http://pyviz.org) shows how to install and work with HoloViews and Datashader together.\n",
    "\n",
    "<style>.container { width:100% !important; }</style>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import datashader as ds\n",
    "import numpy as np\n",
    "import holoviews as hv\n",
    "\n",
    "from holoviews import opts\n",
    "\n",
    "from holoviews.operation.datashader import datashade, shade, dynspread, rasterize\n",
    "from holoviews.operation import decimate\n",
    "\n",
    "hv.extension('bokeh','matplotlib')\n",
    "\n",
    "decimate.max_samples=1000\n",
    "dynspread.max_px=20\n",
    "dynspread.threshold=0.5\n",
    "\n",
    "def random_walk(n, f=5000):\n",
    "    \"\"\"Random walk in a 2D space, smoothed with a filter of length f\"\"\"\n",
    "    xs = np.convolve(np.random.normal(0, 0.1, size=n), np.ones(f)/f).cumsum()\n",
    "    ys = np.convolve(np.random.normal(0, 0.1, size=n), np.ones(f)/f).cumsum()\n",
    "    xs += 0.1*np.sin(0.1*np.array(range(n-1+f))) # add wobble on x axis\n",
    "    xs += np.random.normal(0, 0.005, size=n-1+f) # add measurement noise\n",
    "    ys += np.random.normal(0, 0.005, size=n-1+f)\n",
    "    return np.column_stack([xs, ys])\n",
    "\n",
    "def random_cov():\n",
    "    \"\"\"Random covariance for use in generating 2D Gaussian distributions\"\"\"\n",
    "    A = np.random.randn(2,2)\n",
    "    return np.dot(A, A.T)\n",
    "\n",
    "def time_series(T = 1, N = 100, mu = 0.1, sigma = 0.1, S0 = 20):  \n",
    "    \"\"\"Parameterized noisy time series\"\"\"\n",
    "    dt = float(T)/N\n",
    "    t = np.linspace(0, T, N)\n",
    "    W = np.random.standard_normal(size = N) \n",
    "    W = np.cumsum(W)*np.sqrt(dt) # standard brownian motion\n",
    "    X = (mu-0.5*sigma**2)*t + sigma*W \n",
    "    S = S0*np.exp(X) # geometric brownian motion\n",
    "    return S"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<center><div class=\"alert alert-info\" role=\"alert\">This notebook makes use of dynamic updates, which require a running a live Jupyter or Bokeh server.<br>\n",
    "When viewed statically, the plots will not update fully when you zoom and pan.<br></div></center>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Principles of datashading\n",
    "\n",
    "Because HoloViews elements are fundamentally data containers, not visualizations, you can very quickly declare elements such as ``Points`` or ``Path`` containing datasets that may be as large as the full memory available on your machine (or even larger if using Dask dataframes). So even for very large datasets, you can easily  specify a data structure that you can work with for making selections, sampling, aggregations, and so on. However, as soon as you try to visualize it directly with either the matplotlib or bokeh plotting extensions, the rendering process may be prohibitively expensive.\n",
    "\n",
    "Let's start with a simple example we can visualize as normal:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "np.random.seed(1)\n",
    "points = hv.Points(np.random.multivariate_normal((0,0), [[0.1, 0.1], [0.1, 1.0]], (1000,)),label=\"Points\")\n",
    "paths = hv.Path([random_walk(2000,30)], label=\"Paths\")\n",
    "\n",
    "points + paths"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "These browser-based plots are fully interactive, as you can see if you select the Wheel Zoom or Box Zoom tools and use your scroll wheel or click and drag.  \n",
    "\n",
    "Because all of the data in these plots gets transferred directly into the web browser, the interactive functionality will be available even on a static export of this figure as a web page. Note that even though the visualization above is not computationally expensive, even with just 1000 points as in the scatterplot above, the plot already suffers from [overplotting](https://anaconda.org/jbednar/plotting_pitfalls), with later points obscuring previously plotted points.  \n",
    "\n",
    "With much larger datasets, these issues will quickly make it impossible to see the true structure of the data.  We can easily declare 50X or 1000X larger versions of the same plots above, but if we tried to visualize them they would be nearly unusable even if the browser did not crash:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "np.random.seed(1)\n",
    "points = hv.Points(np.random.multivariate_normal((0,0), [[0.1, 0.1], [0.1, 1.0]], (1000000,)),label=\"Points\")\n",
    "paths = hv.Path([0.15*random_walk(100000) for i in range(10)],label=\"Paths\")\n",
    "\n",
    "#points + paths  ## Danger! Browsers can't handle 1 million points!"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Luckily, HoloViews Elements are just containers for data and associated metadata, not plots, so HoloViews can generate entirely different types of visualizations from the same data structure when appropriate.  For instance, in the plot on the left below you can see the result of applying a `decimate()` operation acting on the `points` object, which will automatically downsample this million-point dataset to at most 1000 points at any time as you zoom in or out:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "decimate(points) + datashade(points) + datashade(paths)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Decimating a plot in this way can be useful, but it discards most of the data, yet still suffers from overplotting. If you have Datashader installed, you can instead use the `datashade()` operation to create a dynamic Datashader-based Bokeh plot. The middle plot above shows the result of using `datashade()` to create a dynamic Datashader-based plot out of an Element with arbitrarily large data.  In the Datashader version, a new image is regenerated automatically on every zoom or pan event, revealing all the data available at that zoom level and avoiding issues with overplotting by dynamically rescaling the colors used.  The same process is used for the line-based data in the Paths plot.\n",
    "\n",
    "These two Datashader-based plots are similar to the native Bokeh plots above, but instead of making a static Bokeh plot that embeds points or line segments directly into the browser, HoloViews sets up a Bokeh plot with dynamic callbacks that render the data as an RGB image using Datashader instead.  The dynamic re-rendering provides an interactive user experience even though the data itself is never provided directly to the browser.  Of course, because the full data is not in the browser, a static export of this page (e.g. on holoviews.org or on anaconda.org) will only show the initially rendered version, and will not update with new images when zooming as it will when there is a live Python process available.\n",
    "\n",
    "Though you can no longer have a completely interactive exported file, with the Datashader version on a live server you can now change the number of data points from 1000000 to 10000000 or more to see how well your machine will handle larger datasets. It will get a bit slower, but if you have enough memory, it should still be very usable, and should never crash your browser as transferring the whole dataset into your browser would.  If you don't have enough memory, you can instead set up a [Dask](http://dask.pydata.org) dataframe as shown in other Datashader examples, which will provide out-of-core and/or distributed processing to handle even the largest datasets.\n",
    "\n",
    "The `datashade()` operation is actually a \"macro\" or shortcut that combines the two main computations done by datashader, namely `shade()` and `rasterize()`:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "rasterize(points).hist() + shade(rasterize(points)) + datashade(points)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In all three of the above plots, `rasterize()` is being called to aggregate the data (a large set of x,y locations) into a rectangular grid, with each grid cell counting up the number of points that fall into it.  In the plot on the left, only `rasterize()` is done, and the resulting numeric array of counts is passed to Bokeh for colormapping.  Bokeh can then use dynamic (client-side, browser-based) operations in JavaScript, allowing users to have dynamic control over  even static HTML plots.  For instance, in this case, users can use the Box Select tool and select a range of the histogram shown, dynamically remapping the colors used in the plot to cover the selected range.\n",
    "\n",
    "The other two plots should be identical.  In both cases, the numerical array output of `rasterize()` is mapped into RGB colors by Datashader itself, in Python (\"server-side\"), which allows special Datashader computations like the histogram-equalization in the above plots and the \"spreading\" discussed below.  The `shade()` and `datashade()` operations accept a `cmap` argument that lets you control the colormap used, which can be selected to match the HoloViews/Bokeh `cmap` option but is strictly independent of it.  See ``hv.help(rasterize)``,  ``hv.help(shade)``, and ``hv.help(datashade)`` for options that can be selected, and the [Datashader web site](http://datashader.org) for all the details. The lower-level `aggregate()` and `regrid()` give more control over how the data is aggregated.\n",
    "\n",
    "Since datashader only sends the data currently in view to the plotting backend, the default behavior is to rescale colormap to the range of the visible data as the zoom level changes. This behavior may not be desirable when working with images; to instead use a fixed colormap range, the `clim` parameter can be passed to the `bokeh` backend via the `opts()` method. Note that this approach works with `rasterize()` where the colormapping is done by the `bokeh` backend. With `datashade()`, the colormapping is done with the `shade()` function which takes a `clims` parameter directly instead of passing additional parameters to the backend via `opts()`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "n = 10_000\n",
    "\n",
    "# Strong signal on top\n",
    "rs = np.random.RandomState(101010)\n",
    "x = rs.pareto(n, n)\n",
    "y = x + rs.standard_normal(n)\n",
    "img1, *_ = np.histogram2d(x, y, bins=60)\n",
    "\n",
    "# Weak signal in the middle\n",
    "x2 = rs.standard_normal(n)\n",
    "y2 = 5 * x + 10 * rs.standard_normal(n)\n",
    "img2, *_ = np.histogram2d(x2, y2, bins=60)\n",
    "\n",
    "img = img1 + img2\n",
    "hv_img = hv.Image(img).opts(active_tools=['wheel_zoom'])\n",
    "auto_scale_grid = rasterize(hv_img).opts(title='Automatic color range rescaling')\n",
    "fixed_scale_grid = rasterize(hv_img).opts(title='Fixed color range', clim=(img.min(), img.max()))\n",
    "auto_scale_grid + fixed_scale_grid; # Output supressed and gif shown below instead"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<img src=\"http://assets.holoviews.org/gifs/guides/user_guide/Large_Data/rasterize_color_range.gif\"></img>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Spreading\n",
    "\n",
    "The Datashader examples above treat points and lines as infinitesimal in width, such that a given point or small bit of line segment appears in at most one pixel. This approach ensures that the overall distribution of the points will be mathematically well founded -- each pixel will scale in value directly by the number of points that fall into it, or by the lines that cross it.\n",
    "\n",
    "However, many monitors are sufficiently high resolution that the resulting point or line can be difficult to see---a single pixel may not actually be visible on its own, and its color may likely be very difficult to make out.  To compensate for this, HoloViews provides access to Datashader's image-based \"spreading\", which makes isolated pixels \"spread\" into adjacent ones for visibility.  There are two varieties of spreading supported:\n",
    "\n",
    "1. ``spread``: fixed spreading of a certain number of pixels, which is useful if you want to be sure how much spreading is done regardless of the properties of the data.\n",
    "2. ``dynspread``: spreads up to a maximum size as long as it does not exceed a specified fraction of adjacency between pixels.  \n",
    "\n",
    "Dynamic spreading is typically more useful, because it adjusts depending on how close the datapoints are to each other on screen. Both types of spreading require Datashader to do the colormapping (applying `shade`), because they operate on RGB pixels, not data arrays.\n",
    "\n",
    "You can compare the results in the two plots below after zooming in:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "datashade(points) + dynspread(datashade(points))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Both plots show the same data, and look identical when zoomed out, but when zoomed in enough you should be able to see the individual data points on the right while the ones on the left are barely visible.  The dynspread parameters typically need some hand tuning, as the only purpose of such spreading is to make things visible on a particular monitor for a particular observer; the underlying mathematical operations in Datashader do not normally need parameters to be adjusted.\n",
    "\n",
    "The same operation works similarly for line segments:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "datashade(paths) + dynspread(datashade(paths))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Multidimensional plots\n",
    "\n",
    "The above plots show two dimensions of data plotted along *x* and *y*, but Datashader operations can be used with additional dimensions as well.  For instance, an extra dimension (here called `k`), can be treated as a category label and used to colorize the points or lines, aggregating the data points separately depending on which category value they have. Compared to a standard overlaid scatterplot that would suffer from overplotting, here the result will be merged mathematically by Datashader, completely avoiding any overplotting issues except any local issues that may arise from spreading when zoomed in:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "np.random.seed(3)\n",
    "kdims=['d1','d2']\n",
    "num_ks=8\n",
    "\n",
    "def rand_gauss2d(value=0, n=100000):\n",
    "    \"\"\"Return a randomly shaped 2D Gaussian distribution with an associated numeric value\"\"\"\n",
    "    g = 100*np.random.multivariate_normal(np.random.randn(2), random_cov(), (n,))\n",
    "    return np.hstack((g,value*np.ones((g.shape[0],1))))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "gaussians = {i: hv.Points(rand_gauss2d(i), kdims, \"i\") for i in range(num_ks)}\n",
    "c = dynspread(datashade(hv.NdOverlay(gaussians, kdims='k'), aggregator=ds.by('k', ds.count())))\n",
    "m = dynspread(datashade(hv.NdOverlay(gaussians, kdims='k'), aggregator=ds.by('k', ds.mean(\"i\"))))\n",
    "\n",
    "(c + m).opts(opts.RGB(width=400))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Above you can see that (as of Datashader 0.11) categorical aggregates can take any reduction function, either `count`ing the datapoints (left) or reporting some other statistic (e.g. the mean value of a column, right).\n",
    "\n",
    "Categorical aggregates are one way to allow separate lines or other shapes to be visually distinctive from one another while avoiding obscuring data due to overplotting:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "lines = {i: hv.Curve(time_series(N=10000, S0=200+np.random.rand())) for i in range(num_ks)}\n",
    "linespread = dynspread(datashade(hv.NdOverlay(lines, kdims='k'), aggregator=ds.by('k', ds.count())))\n",
    "linespread.opts(opts.RGB(width=400))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note that Bokeh only ever sees an image come out of `hd.datashade`, not any of the actual data. As a result, providing legends and keys has to be done separately, though we are working to make this process more seamless.  For now, you can show a legend by adding a suitable collection of \"fake\" labeled points (size zero and thus invisible):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# definition copied here to ensure independent pan/zoom state for each dynamic plot\n",
    "gaussspread2 = dynspread(datashade(hv.NdOverlay(gaussians, kdims=['k']), aggregator=ds.by('k', ds.count())))\n",
    "\n",
    "from datashader.colors import Sets1to3 # default datashade() and shade() color cycle\n",
    "color_key = list(enumerate(Sets1to3[0:num_ks]))\n",
    "color_points = hv.NdOverlay({k: hv.Points([0,0], label=str(k)).opts(color=v, size=0) for k, v in color_key})\n",
    "\n",
    "(color_points * gaussspread2).opts(width=600)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Here the dummy points are at 0,0 for this dataset, but would need to be at another suitable value for data that is in a different range."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Working with time series"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "HoloViews also makes it possible to datashade large timeseries using the ``datashade`` and ``rasterize`` operations:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "\n",
    "dates = pd.date_range(start=\"2014-01-01\", end=\"2016-01-01\", freq='1D') # or '1min'\n",
    "curve = hv.Curve((dates, time_series(N=len(dates), sigma = 1)))\n",
    "datashade(curve, cmap=[\"blue\"], width=800).opts(width=800)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "HoloViews also supplies some operations that are useful in combination with Datashader timeseries.  For instance, you can compute a rolling mean of the results and then show a subset of outlier points, which will then support hover, selection, and other interactive Bokeh features:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from holoviews.operation.timeseries import rolling, rolling_outlier_std\n",
    "smoothed = rolling(curve, rolling_window=50)\n",
    "outliers = rolling_outlier_std(curve, rolling_window=50, sigma=2)\n",
    "\n",
    "ds_curve = datashade(curve, cmap=[\"blue\"])\n",
    "spread = dynspread(datashade(smoothed, cmap=[\"red\"], width=800),max_px=1) \n",
    "\n",
    "(ds_curve * spread * outliers).opts(\n",
    "    opts.Scatter(line_color=\"black\", fill_color=\"red\", size=10, tools=['hover', 'box_select'], width=800))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note that the above plot will look blocky in a static export (such as on anaconda.org), because the exported version is generated without taking the size of the actual plot (using default height and width for Datashader) into account, whereas the live notebook automatically regenerates the plot to match the visible area on the page. The result of all these operations can be laid out, overlaid, selected, and sampled just like any other HoloViews element, letting you work naturally with even very large datasets.\n",
    "\n",
    "\n",
    "# Hover info\n",
    "\n",
    "As you can see in the examples above, converting the data to an image using Datashader makes it feasible to work with even very large datasets interactively.  One unfortunate side effect is that the original datapoints and line segments can no longer be used to support \"tooltips\" or \"hover\" information directly for RGB images generated with `datashade`; that data simply is not present at the browser level, and so the browser cannot unambiguously report information about any specific datapoint. \n",
    "\n",
    "If you do need hover information, there are two good options available:\n",
    "\n",
    "1) Use the ``rasterize`` operation without `shade`, which will let the plotting code handle the conversion to colors while still having the actual aggregated data to support hovering\n",
    "\n",
    "2) Overlay a separate layer as a ``QuadMesh`` or ``Image`` containing the hover information"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from holoviews.streams import RangeXY\n",
    "\n",
    "rasterized = rasterize(points, width=400, height=400)\n",
    "\n",
    "fixed_hover = (datashade(points, width=400, height=400) *  \n",
    "               hv.QuadMesh(rasterize(points, width=10, height=10, dynamic=False)))\n",
    "\n",
    "dynamic_hover = (datashade(points, width=400, height=400) * \n",
    "                 rasterize(points, width=10, height=10, streams=[RangeXY]).apply(hv.QuadMesh))\n",
    "\n",
    "(rasterized + fixed_hover + dynamic_hover).opts(\n",
    "    opts.QuadMesh(tools=['hover'], alpha=0, hover_alpha=0.2), \n",
    "    opts.Image(tools=['hover']))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In the above examples, the plot on the left provides hover information directly on the aggregated ``Image``. The middle plot displays hover information as a ``QuadMesh`` at a fixed spatial scale, while the one on the right reports on an area that scales with the zoom level so that arbitrarily small regions of data space can be examined, which is generally more useful (but requires a live Python server)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Element types supported for Datashading\n",
    "\n",
    "Fundamentally, what datashader does is to rasterize data, i.e., render a representation of it into a regularly gridded rectangular portion of a two-dimensional plane.  Datashader natively supports six basic types of rasterization:\n",
    "\n",
    "- **points**: zero-dimensional objects aggregated by position alone, each point covering zero area in the plane and thus falling into exactly one grid cell of the resulting array (if the point is within the bounds being aggregated).\n",
    "- **line**: polyline/multiline objects (connected series of line segments), with each segment having a fixed length but zero width and crossing each grid cell at most once.\n",
    "- **area**: a region to fill either between the supplied y-values and the origin or between two supplied lines\n",
    "- **trimesh**: irregularly spaced triangular grid, with each triangle covering a portion of the 2D plane and thus potentially crossing multiple grid cells (thus requiring interpolation/upsampling). Depending on the zoom level, a single pixel can also include multiple triangles, which then becomes similar to the `points` case (requiring aggregation/downsampling of all triangles covered by the pixel).\n",
    "- **raster**: an axis-aligned regularly gridded two-dimensional subregion of the plane, with each grid cell in the source data covering more than one grid cell in the output grid (requiring interpolation/upsampling), or with each grid cell in the output grid including contributions from more than one grid cell in the input grid (requiring aggregation/downsampling).\n",
    "- **quadmesh**: a recti-linear or curvi-linear mesh (like a raster, but allowing nonuniform spacing and coordinate mapping) where each quad can cover one or more cells in the output (requiring upsampling, currently only as nearest neighbor), or with each output grid cell including contributions from more than one input grid cell (requiring aggregation/downsampling).\n",
    "- **polygons**:  arbitrary filled shapes in 2D space (bounded by a piecewise linear set of segments), optionally punctuated by similarly bounded internal holes\n",
    "\n",
    "Datashader focuses on implementing those four cases very efficiently, and HoloViews in turn can use them to render a very large range of specific types of data: \n",
    "\n",
    "### Supported Elements\n",
    "\n",
    "- **points**: [`hv.Nodes`](../reference/elements/bokeh/Graph.ipynb), [`hv.Points`](../reference/elements/bokeh/Points.ipynb), [`hv.Scatter`](../reference/elements/bokeh/Scatter.ipynb)\n",
    "- **line**: [`hv.Contours`](../reference/elements/bokeh/Contours.ipynb), [`hv.Curve`](../reference/elements/bokeh/Curve.ipynb), [`hv.Path`](../reference/elements/bokeh/Path.ipynb), [`hv.Graph`](../reference/elements/bokeh/Graph.ipynb), [`hv.EdgePaths`](../reference/elements/bokeh/Graph.ipynb), [`hv.Spikes`](../reference/elements/bokeh/Spikes.ipynb), [`hv.Segments`](../reference/elements/bokeh/Segments.ipynb)\n",
    "- **area**: [`hv.Area`](../reference/elements/bokeh/Area.ipynb), [`hv.Spread`](../reference/elements/bokeh/Spread.ipynb)\n",
    "- **raster**: [`hv.Image`](../reference/elements/bokeh/Image.ipynb), [`hv.HSV`](../reference/elements/bokeh/HSV.ipynb), [`hv.RGB`](../reference/elements/bokeh/RGB.ipynb)\n",
    "- **trimesh**: [`hv.TriMesh`](../reference/elements/bokeh/TriMesh.ipynb)\n",
    "- **quadmesh**: [`hv.QuadMesh`](../reference/elements/bokeh/QuadMesh.ipynb)\n",
    "- **polygons**: [`hv.Polygons`](../reference/elements/bokeh/Polygons.ipynb)\n",
    "\n",
    "Other HoloViews elements *could* be supported, but do not currently have a useful datashaded representation:\n",
    "\n",
    "### Elements not yet supported\n",
    "\n",
    "- **line**: [`hv.Spline`](../reference/elements/bokeh/Spline.ipynb), [`hv.VectorField`](../reference/elements/bokeh/VectorField.ipynb)\n",
    "- **raster**: [`hv.HeatMap`](../reference/elements/bokeh/HeatMap.ipynb), [`hv.Raster`](../reference/elements/bokeh/Raster.ipynb)\n",
    "\n",
    "There are also other Elements that are not expected to be useful with datashader because they are isolated annotations, are already summaries or aggregations of other data, have graphical representations that are only meaningful at a certain size, or are text based:\n",
    "\n",
    "### Not useful to support\n",
    "\n",
    "- datashadable annotations: [`hv.Arrow`](../reference/elements/bokeh/Arrow.ipynb), [`hv.Bounds`](../reference/elements/bokeh/Bounds.ipynb), [`hv.Box`](../reference/elements/bokeh/Box.ipynb), [`hv.Ellipse`](../reference/elements/bokeh/Ellipse.ipynb) (actually do work with datashade currently, but not officially supported because they are not vectorized and thus unlikely to have enough items to be worth datashading)\n",
    "- other annotations: [`hv.Arrow`](../reference/elements/bokeh/Arrow.ipynb), [`hv.HLine`](../reference/elements/bokeh/HLine.ipynb), [`hv.VLine`](../reference/elements/bokeh/VLine.ipynb), [`hv.Text`](../reference/elements/bokeh/Text.ipynb)\n",
    "- kdes: [`hv.Distribution`](../reference/elements/bokeh/Distribution.ipynb), [`hv.Bivariate`](../reference/elements/bokeh/Bivariate.ipynb)\n",
    "- categorical/symbolic: [`hv.BoxWhisker`](../reference/elements/bokeh/BoxWhisker.ipynb), [`hv.Bars`](../reference/elements/bokeh/Bars.ipynb), [`hv.ErrorBars`](../reference/elements/bokeh/ErrorBars.ipynb)\n",
    "- tables: [`hv.Table`](../reference/elements/bokeh/Table.ipynb), [`hv.ItemTable`](../reference/elements/bokeh/ItemTable.ipynb)\n",
    "\n",
    "Examples of each supported Element type:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "hv.output(backend='matplotlib')\n",
    "\n",
    "from bokeh.sampledata.us_counties import data as counties\n",
    "\n",
    "opts.defaults(\n",
    "    opts.Image(aspect=1, axiswise=True, xaxis='bare', yaxis='bare'),\n",
    "    opts.RGB(aspect=1, axiswise=True, xaxis='bare', yaxis='bare'),\n",
    "    opts.HSV(aspect=1, axiswise=True, xaxis='bare', yaxis='bare'),\n",
    "    opts.Layout(vspace=0.1, hspace=0.1, sublabel_format=\"\", fig_size=80))\n",
    "\n",
    "np.random.seed(12)\n",
    "N=100\n",
    "pts = [(10*i/N, np.sin(10*i/N)) for i in range(N)]\n",
    "\n",
    "x = y = np.linspace(0, 5, int(np.sqrt(N)))\n",
    "xs,ys = np.meshgrid(x,y)\n",
    "z = np.sin(xs)*np.cos(ys)\n",
    "\n",
    "r = 0.5*np.sin(0.1*xs**2+0.05*ys**2)+0.5\n",
    "g = 0.5*np.sin(0.02*xs**2+0.2*ys**2)+0.5\n",
    "b = 0.5*np.sin(0.02*xs**2+0.02*ys**2)+0.5\n",
    "\n",
    "n=20\n",
    "coords = np.linspace(-1.5,1.5,n)\n",
    "X,Y = np.meshgrid(coords, coords);\n",
    "Qx = np.cos(Y) - np.cos(X)\n",
    "Qy = np.sin(Y) + np.sin(X)\n",
    "Z = np.sqrt(X**2 + Y**2)\n",
    "\n",
    "s = np.random.randn(100).cumsum()\n",
    "e = s + np.random.randn(100)\n",
    "\n",
    "opts2 = dict(filled=True, edge_color='z')\n",
    "tri = hv.TriMesh.from_vertices(hv.Points(np.random.randn(N,3), vdims='z')).opts(**opts2)\n",
    "(tri + tri.edgepaths + datashade(tri, aggregator=ds.mean('z')) + datashade(tri.edgepaths)).cols(2)\n",
    "\n",
    "shadeable  = [elemtype(pts) for elemtype in [hv.Curve, hv.Scatter]]\n",
    "shadeable += [hv.Path(counties[(1, 1)], ['lons', 'lats']), hv.Points(counties[(1, 1)], ['lons', 'lats'])]\n",
    "shadeable += [hv.Spikes(np.random.randn(10000))]\n",
    "shadeable += [hv.Segments((np.arange(100), s, np.arange(100), e))]\n",
    "shadeable += [hv.Area(np.random.randn(10000).cumsum())]\n",
    "shadeable += [hv.Spread((np.arange(10000), np.random.randn(10000).cumsum(), np.random.randn(10000)*10))]\n",
    "shadeable += [hv.Image((x,y,z))]\n",
    "shadeable += [hv.QuadMesh((Qx,Qy,Z))]\n",
    "shadeable += [hv.Graph(((np.zeros(N), np.arange(N)),))]\n",
    "shadeable += [tri.edgepaths]\n",
    "shadeable += [tri]\n",
    "shadeable += [hv.Polygons([county for county in counties.values() if county['state'] == 'tx'], ['lons', 'lats'], ['name'])]\n",
    "shadeable += [hv.operation.contours(hv.Image((x,y,z)), levels=10)]\n",
    "\n",
    "rasterizable = [hv.RGB(np.dstack([r,g,b])), hv.HSV(np.dstack([g,b,r]))]\n",
    "\n",
    "ds_opts = {\n",
    "    hv.Path: dict(aggregator='any'),\n",
    "    hv.Points: dict(aggregator='any'),\n",
    "    hv.Polygons: dict(aggregator=ds.count_cat('name'), color_key=hv.plotting.util.process_cmap('glasbey')),\n",
    "    hv.Segments: dict(aggregator='any')\n",
    "}\n",
    "\n",
    "hv.Layout([dynspread(datashade(e.relabel(e.__class__.name), **ds_opts.get(e.__class__, {}))) for e in shadeable] + \n",
    "          [          rasterize(e.relabel(e.__class__.name))  for e in rasterizable]).opts(shared_axes=False).cols(6)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Here we called `datashade()` on each Element type, letting Datashader do the full process of rasterization and shading, except that for `RGB` and `HSV` we only called `rasterize()` or else the results would have been converted into a monochrome image.\n",
    "\n",
    "For comparison, you can see the corresponding non-datashaded plots (as long as you leave N lower than 10000 unless you have a long time to wait!):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "el_opts = dict(aspect=1, axiswise=True, xaxis='bare', yaxis='bare')\n",
    "hv.Layout([e.relabel(e.__class__.name).opts(**el_opts) for e in shadeable + rasterizable]).cols(6)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "These two examples use Matplotlib, but if they were switched to Bokeh and you had a live server, they would support dynamic re-rendering on zoom and pan so that you could explore the full range of data available (e.g. even very large raster images, networks, paths, point clouds, or meshes).\n",
    "\n",
    "\n",
    "# Container types supported for datashading\n",
    "\n",
    "In the above examples `datashade()` was called directly on each Element, but it can also be called on Containers, in which case each Element in the Container will be datashaded separately (for all Container types other than a Layout):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "hv.output(dpi=80, size=100)\n",
    "\n",
    "curves = {'+':hv.Curve(pts), '-':hv.Curve([(x, -1.0*y) for x, y in pts])}\n",
    "\n",
    "supported = [hv.HoloMap(curves,'sign'), hv.Overlay(list(curves.values())), hv.NdOverlay(curves), hv.GridSpace(hv.NdOverlay(curves))]\n",
    "hv.Layout([datashade(e.relabel(e.__class__.name)) for e in supported]).cols(4)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "dynspread(datashade(hv.NdLayout(curves,'sign')))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "hv.output(backend='bokeh')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Optimizing performance\n",
    "\n",
    "Datashader and HoloViews have different design principles that are worth keeping in mind when using them in combination, if you want to ensure good overall performance. By design, Datashader supports only a small number of operations and datatypes, focusing only on what can be implemented very efficiently.  HoloViews instead focuses on supporting the typical workflows of Python users, recognizing that the most computationally efficient choice is only going to be faster overall if it also minimizes the time users have to spend getting things working. \n",
    "\n",
    "HoloViews thus helps you get something working quickly, but once it is working and you realize that you need to do this often or that it comes up against the limits of your computing hardware, you can consider whether you can get much better performance by considering the following issues and suggestions.\n",
    "\n",
    "### Use a Datashader-supported data structure\n",
    "\n",
    "HoloViews helpfully tries to convert whatever data you have provided into what Datashader supports, which is good for optimizing your time to an initial solution, but will not always be the fastest approach computationally.  If you ensure that you store your data in a format that Datashader understands already, HoloViews can simply pass it down to Datashader without copying or transforming it:\n",
    "\n",
    "1. For point, line, and trimesh data, Datashader supports Dask and Pandas dataframes, and so those two data sources will be fastest. Of those two, Dask Dataframes will usually be somewhat faster and also able to make use of distributed computational resources and out-of-core processing.\n",
    "2. For rasters and quadmeshes, Datashader supports xarray objects natively, and so if your data is provided as an xarray already, plotting will be faster.\n",
    "3. For polygons Datashader supports [spatialpandas](https://github.com/holoviz/spatialpandas) DataFrames.\n",
    "\n",
    "See the [Datashader docs](http://datashader.org) for examples of dealing with even quite large datasets (in the billions of points) on commodity hardware, including many HoloViews-based examples.\n",
    "\n",
    "### Cache initial processing with `precompute=True`\n",
    "\n",
    "In the typical case of having datasets much larger than the plot resolution, HoloViews Datashader-based operations that work on the full dataset (`rasterize`, `aggregate`,`regrid`) are computationally expensive; the others are not (`shade`, `spread`, `dynspread`, etc.) \n",
    "\n",
    "The expensive operations are all of type `ResamplingOperation`, which has a parameter `precompute` (see `hv.help(hv.operation.datashader.rasterize)`, etc.)  Precompute can be used to get faster performance in interactive usage by caching the last set of data used in plotting (*after* any transformations needed) and reusing it when it is requested again. This is particularly useful when your data is not in one of the supported data formats already and needs to be converted. `precompute` is False by default, because it requires using memory to store the cached data, but if you have enough memory, you can enable it so that repeated interactions (such as zooming and panning) will be much faster than the first one.  In practice, most Datashader-plots don't need to do extensive precomputing, but enabling it for TriMesh and Polygon plots can greatly speed up interactive usage.\n",
    "\n",
    "### Use GPU support\n",
    "\n",
    "Many elements now also support aggregation directly on a GPU-based datastructure such as a [cuDF DataFrame](https://github.com/rapidsai/cudf) or an Xarray DataArray backed by a [cupy](https://github.com/cupy/cupy) array. These data structures can be passed directly to the appropriate HoloViews elements just as you would use a Pandas or other Xarray object. For instance, a cuDF can be used on elements like `hv.Points` and `hv.Curve`, while a cupy-backed DataArray raster or quadmesh can be passed to `hv.QuadMesh` elements. When used with Datashader, the GPU implementation can result in 10-100x speedups, as well as avoiding having to transfer the data out of the GPU for plotting (sending only the final rendered plot out of the GPU's memory). To see which HoloViews elements are supported, see the [datashader performance guide](https://datashader.org/user_guide/Performance.html). As of the Datashader 0.11 release, all point, line, area, and quadmesh aggregations are supported when using a GPU backed datastructure, including raster objects like `hv.Image` if first converted to `hv.Quadmesh`.\n",
    "\n",
    "### Project data only once\n",
    "\n",
    "If you are working with geographic data using [GeoViews](http://geoviews.org) that needs to be projected before display and/or before datashading, GeoViews will have to do this every time you update a plot, which can drown out the performance improvement you get by using Datashader.  GeoViews allows you to project the entire dataset at once using `gv.operation.project`, and once you do this you should be able to use Datashader at full speed.\n",
    "\n",
    "If you follow these suggestions, the combination of HoloViews and Datashader will allow you to work uniformly with data covering a huge range of sizes. Per session or per plot, you can trade off the ability to export user-manipulable plots against file size and browser compatibility, and allowing you to render even the largest dataset faithfully. HoloViews makes the full power of Datashader available in just a few lines of code, giving you a natural way to work with your data regardless of its size."
   ]
  }
 ],
 "metadata": {
  "language_info": {
   "name": "python",
   "pygments_lexer": "ipython3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}