By Philippe Saadé and Olivier Cruchant

Introduction

Apache MXNet is a modern deep learning framework featuring both a pleasant developer experience and high-performance training and inference. In order to get the best possible performance while exposing a friendly interface, the MXNet community re-developed from scratch a suite of optimized tensor algebra kernels, historically accessible via the Python MXNet NDArray library. The initial MXNet NDArray API is very similar to NumPy, yet it features major internal differences enabling high-performance deep learning: (1) it can run on GPU, (2) it supports auto-differentiation, (3) it runs asynchronously, executing only the code that need to run and in an optimized order. In order to make the powerful MXNet NDArray API fully accessible to NumPy developers, the MXNet community released in 2019 the mxnet.numpy library, which implements the NumPy API (Announcement on MXNet Medium). In this blog post, we illustrate the strength of MXNet NDArray library for fast and compact algebra over images. We first highlight couple NDArray concepts and then assemble them into a function, that we provide both as an mxnet.ndarray implementation and an mxnet.numpy implementation. We use as an illustration a fictional use-case of anomalous image detection via simple pixel difference. This can be a reasonable baseline when looking for anomalous areas over a batch of near-identical images with same viewpoint and luminosity, for example the images below:

Notable concepts

We will use the following libraries:

import itertools as itr
import os
import mxnet as mx
from mxnet import image as mxim
from mxnet import init, gluon, nd
import matplotlib.cm as cm
from matplotlib import pyplot as plt
from matplotlib import image as mpim

MXNet features the image library that we can use to read images:

picture = mxim.imread('orig1.JPG')

# .asnumpy() converts MXNet NDArrays to Numpy array
plt.imshow(picture.asnumpy())

MXNet NDArray leverages the mxnet context concept, that enables to control the hardware context - CPU or GPU of objects. For example, we can send our previously loaded image to the GPU (if on a GPU-equipped instance) with the following snippet:

pic_on_gpu = picture.copyto(mx.gpu())

Note that we can also use picture.as_in_context(mx.gpu()) to copy the picture to the GPU ; the advantage of the latter being that it does the actual copy only if the image is not in the GPU already.

Similarly to Numpy, MXNet NDArray is a great platform to apply algebraic transformations to batches of records. Writing transformations applied to batch of images is furthermore relevant on GPUs which excel at batch processing. Very similar to their Numpy counterparts, we use in this demo ndarray.concat(), ndarray.mean(), ndarray.abs(), ndarray.max() to manipulate and analyze a batch of images.

MXNet being primarily designed for deep learning, it features abundant neural network primitives. Its imperative Python front-end Gluon provides numerous model layers, and in this blog - even though we do not use deep learning nor even machine learning - we borrow from Gluon its 2D-convolution (mxnet.gluon.nn.Conv2D) with all coefficients initialized at 1 to apply a smoothing effect and isolate contiguous areas of anomalous pixels.

Wrapping it all in a compact function

The inline code block below proposes an example function to run the above-mentioned analysis. It first reads a batch of images, then subtracts from all of them the cross-batch mean image, applies a hard-coded convolution to each of them and return images that have high-value pixels after the convolution. Note that this code features couple hard-coded constants that would be worth tuning over representative data. This sample is provided as an illustration of the NDArray capacities and may deserve further refinement and testing before being used in the real-world.

Implementation with mxnet.ndarray

def mxnd_find_anomalies(folder, gpu=False, save_viz=True, threshold=300):
    """
    Parameters
        ----------
        folder : str
            local directory with the batch of images to score
        gpu : bool, optional
            whether to use GPU or CPU
        save_viz : bool, optional
            saves the pictures of deltas and convolutions
        threshold : int, optional
            conv threshold to be anomalous. Should be tuned!           
    """
    
    ctx = mx.gpu() if gpu else mx.cpu()  # Set context
    
    # read images in a 4D-NDArray
    pics = os.listdir(folder)
    ims = [mxim.imread(folder+'/'+pic).expand_dims(0).as_in_context(ctx) for pic in pics]
    ims = nd.concat(*ims, dim=0)
        
    # Compute the average image
    avg = nd.mean(ims.astype('float32'), axis=0)
    
    # Remove the mean from every image, then average-pool on color
    deltas = nd.mean(ims.astype('float32')-avg, axis=3)
    
    # Apply a batched convolution to all the deltas, then apply abs
    conv = gluon.nn.Conv2D(1, kernel_size=5, use_bias=False)
    conv.initialize(init.Constant(1))
    conv.collect_params().reset_ctx(ctx)
    
    # we expand dims for conv, which expect 3D pics (with channels)
    conv_delta = nd.abs(conv(deltas.expand_dims(1)))
    
    # Max-pool over the image
    top_deltas = nd.max(conv_delta, axis=(1,2,3))
    
    # Return pics that have a convolution > threshold
    anomalies = list(itr.compress(pics, top_deltas > threshold))
    an_indexes = list(itr.compress(range(len(pics)), top_deltas > threshold))
    
    if save_viz:  # optionally, save images for interpretability
        for a, i in zip(anomalies, an_indexes):
            mpim.imsave('pixmap-'+a, deltas[i].asnumpy(), format='png')
            mpim.imsave('convmap-'+a, conv_delta[i][0].asnumpy(), format='png')
                
    return list(anomalies)

Implementation with mxnet.numpy

In the below code block, we implement the same logic as above, while replacing the mxnet.ndarray API by the mxnet.numpy API. Since mxnet.numpy implements the NumPy API, we need to replace couple methods by their NumPy counterpart, notably:

• Expanding array dimensions: array.expand_dims(0) above becomes np.expand_dims(array, 0)
• Concatenation: nd.concat(*ims, dim=0) becomes np.concatenate(ims, axis=0)

We further need to import mxnet.numpy and notify MXNet backend of our use of NumPy semantics, with the extra import and configuration below:

from mxnet import np, npx
npx.set_np()

And here is the final implementation of our baseline anomaly detection function using the MXNet NumPy API:

def mxnp_find_anomalies(folder, gpu=False, save_viz=True, threshold=300):
    """
    Parameters
        ----------
        folder : str
            local directory with the batch of images to score
        gpu : bool, optional
            whether to use GPU or CPU
        save_viz : bool, optional
            saves the pictures of deltas and convolutions
        threshold : int, optional
            conv threshold to be anomalous. Should be tuned!
    """
    
    ctx = mx.gpu() if gpu else mx.cpu()  # Set context
    
    # read images in a 4D-NDArray
    pics = os.listdir(folder)
    ims = [np.expand_dims(mxim.imread(folder+'/'+pic), 0).as_in_context(ctx) for pic in pics]
    ims = np.concatenate(ims, axis=0)
        
    # Compute the average image
    avg = np.mean(ims.astype('float32'), axis=0)
    
    # Remove the mean from every image, then average-pool on color
    deltas = np.mean(ims.astype('float32')-avg, axis=3)
    
    # Apply a batched convolution to all the deltas, then apply abs
    conv = gluon.nn.Conv2D(1, kernel_size=5, use_bias=False)
    conv.initialize(init.Constant(1))
    conv.collect_params().reset_ctx(ctx)
    
    # we expand dims for conv, which expect 3D pics (with channels)
    conv_delta = np.abs(conv(np.expand_dims(deltas, 1)))
    
    # Max-pool over the image
    top_deltas = np.max(conv_delta, axis=(1,2,3))
    
    # Return pics that have a convolution > threshold
    anomalies = list(itr.compress(pics, top_deltas > threshold))
    an_indexes = list(itr.compress(range(len(pics)), top_deltas > threshold))
    
    if save_viz:  # optionally, save images for interpretability
        for a, i in zip(anomalies, an_indexes):
            mpim.imsave('pixmap-'+a, deltas[i].asnumpy(), format='png')
            mpim.imsave('convmap-'+a, conv_delta[i][0].asnumpy(), format='png')
                
    return list(anomalies)

Here are 2 anomalies the function found over our synthetic data:

Both functions have comparable latencies, which makes sense since they use the same back-end kernels. Here are execution latencies excluding visualizations creation, over a batch of 12 1000px*750px images with MXNet 1.7, on an Amazon SageMaker GPU-equipped ml.p3.2xlarge Notebook (average of 10 runs)

• The mxnet.ndarray code runs in 780ms and the mxnet.numpy code in 750ms on the 8 vCPU (Intel Xeon E5)
• The mxnet.ndarray code runs in 100ms and the mxnet.numpy code in 75ms on the GPU (NVIDIA Tesla V100 GPU), close to 90% faster vs CPU

Conclusion

Even though MXNet is a deep learning framework, its primitives can be creatively used beyond neural network development. In this post, we showed how to perform simple, GPU-accelerated algebra over images using its NDArray library, that can be used by NumPy developers without extra learning curve via the mxnet.numpy library. Do not hesitate to take a look at it, contribute and engage with the community!

Using MXNet NDArray for fast GPU algebra on images was originally published in Apache MXNet on Medium, where people are continuing the conversation by highlighting and responding to this story.

It is exciting to see the fast development of commercial services and open-source tooling supporting inference optimization. While the sheer volume of innovation in the space can sometimes feel overwhelming, a simple mental model to understand those optimizations techniques and tricks is to cluster them along three orthogonal dimensions: (1) hardware, (2) software and (3) scientific logic of the deployed algorithm. Eventually, we will also discuss architectural consideration that can contribute to further speed and cost optimizations such as caching, serverless deployment and model locality.

1. Hardware optimizations

A usually straightforward technique to reduce model latency is to use an appropriate hardware back-end. You can turn those bolts in two directions:

Horizontal scaling — There is usually a bidirectional causal relationship between server load and response latency, where latency gets degraded as server load increases, and server load increases if the model takes longer to run. In that regard, to keep latency predictably low, increasing the number of server nodes is a first technique that can be done manually but also programmatically. For example, the Autoscaling feature of Amazon SageMaker endpoints enables you to configure scale-out and scale-in policies based on a metric of your choice, typically using the pre-defined Amazon Cloudwatch metric InvocationsPerInstance. Determining appropriate thresholds for the scaling policy can be done via load testing, a best practice that can be part of the model deployment cycle. Note that some forms of load testing qualify for the Amazon EC2 Testing Policy. An AWS ML Blog post presents the details of SageMaker Autoscaling and end point load testing.

Vertical scaling & hardware specialization - Another axis of hardware back-end improvement is to right-size the endpoint instances and choose appropriate compute host type and AI accelerator, if relevant to the model. The concept of AI accelerator usually refers to the class of microprocessors designed for AI applications. NVIDIA GPUs are a common option for AI acceleration, most notably the recent Volta and Turing architectures that feature the Tensor Cores, which are specialized cells for mixed-precision matrix multiply-accumulate are found in the AWS P3 and G4 families of Amazon EC2 instances respectively. Notably, other AI accelerators are being made available, such as Amazon Elastic Inference (EIA) and AWS Inferentia. Every model, ML framework and serving stack has its sweet spot when it comes to choosing the right mix of RAM, CPU, accelerator memory and accelerator compute capacity. Appropriate load testing will reveal which instance is more appropriate for your use-case. The use of AI accelerators, such as GPU or Amazon Elastic Inference, usually accelerates the inference.

An interesting comparison was recently made by a customer of AWS, Curalate, in the AWS ML Blog, comparing inference of various vision models in AWS Lambda back-end (serverless, CPU) vs Amazon EC2 with Amazon Elastic Inference acceleration. Researchers from Curalate were successful in running Apache MXNet vision classifiers inference in Amazon Lambda, and found for various models thresholds of traffic below which it AWS Lambda deployment would be cheaper than when using always-on EC2. For example, they measured that in their setup, a ResNet 152 deployment would be cheaper than an EC2 + Amazon Elastic Inference-based deployment for traffics up to 7.5MM calls/month. They also mention that leveraging Amazon Elastic Inference led to cost savings on hosting between 35% and 65% compared to regular Amazon EC2 GPU instances.

2. Software optimizations

Another axis of improvement, separated from the hardware and algorithm considerations, is to work on the software stack, itself composed of two sub-parts — server and compilation.

The model server is tasked with loading one or several models into memory and exposing them through a web service. Different servers will handle load at different levels of performance and provide various levels of convenience. A common trick to maximize hosting hardware utilization is to leverage batching. For example, in some situations instead of calling your server 100 times, with 100 individual requests, you can call your server 1 time with a batch of 100 requests. Some ML problems such as ad scoring or recommendation naturally lends themselves to batching, as you can for example score relevancy of all products in a catalog in one model call. Specialized deep learning servers such as TensorFlow Serving, TensorRT Inference Server or Multi Model Server (MMS) are actually able to accumulate incoming requests into batches to improve GPU parallelism, a practice usually called as micro-batching. While micro-batching can improve inference economics, it usually comes at the expense of latency, since incoming requests will wait a varying amount of time while the micro-batch is assembled. This logic is handled by the server, and servers with micro-batch enabled usually come with controls to tune maximum acceptable latency and micro-batch size. While MMS is not specific to Apache MXNet, it was built with strong native integration with Apache MXNet and features a number of functionalities specific to ML model serving, such as logging, model pre-loading, micro-batching and multi-model serving.

Another axis of software-based improvement is model compilation, a post-training step that consists of creating a model artifact and a model execution runtime finely adapted to the underlying supporting hardware. This concept leverages the fact that when you know where you’re going to deploy (e.g., NVIDIA, Intel, ARM), you have an insider edge and you can refine the model graph and inference runtime to better leverage available low-level specifics. This can reduce memory consumption and latency by double-digit percentage, and is an active area of ML research as illustrated by the papers mentioned in the rest of this paragraph. A popular, flexible and high-performance option for compilation is Tensor Virtual Machine (TVM, Chen et al, paper, documentation), which demonstrated compilation from a large number of frameworks (MXNet, Keras, Caffe2, CoreML, TensorFlow,…). Other notable compilers are for example treelite (paper, doc) for decision tree compilation, and the Glow compilation engine for PyTorch. Interestingly, both treelite and TVM are dependencies of the AWS compiler SageMaker Neo, which provides a managed experience of compilation.

One benefit of using Amazon SageMaker Neo is that it is compatible with the most popular frameworks (MXNet, TensorFlow, XGBoost, PyTorch) and hardware back-ends (Intel, NVIDIA, ARM and others), and that it chooses on your behalf the underlying compiler tricks most appropriate for your combination of hardware target, framework and model.
Model quantization refers to reducing the number of bits representing model weights. It is usually available as part of a compilation stack or as a model framework component (examples of quantization functions in TVM, TensorFlow Lite, Apache MXNet or NVIDIA TensorRT). Quantization aims to reduce model storage footprint and computation time, and is most commonly a post-training step. Another benefit of quantization is that it can leverage specific low-precision optimizations, such as the Vector Neural Network Instruction (VNNI) on modern Intel Xeon Scalable CPUs. Quantization is most commonly ran post-training, with an optional calibration step during which the model sees sample data to properly quantify the ranges to be quantized. Quantization may also be run in-training, a concept named quantization-aware training that involves directly training a quantized representation. At the time of this writing, this seems to be still at active research stage and less a popular option than post-training quantization. Notably, in summer 2019 the TensorFlow group encouraged in a blog post the use of post-training quantization over quantization-aware training.

3. Algorithm optimization
On a given hardware and software platform, model inference latency and hardware footprint will significantly vary with the algorithmic logic of the model. If speed, throughput or costs are your priorities over accuracy, and if you already exploited all possible tricks along hardware and software axis above, consider using models adapted to those goals. For example:

• Models with very little computation, such as logistic regression for classification over structured data
• In the computer vision classification space, the popular MobileNet architecture (Howard et al. https://arxiv.org/abs/1704.04861) enables significant computation and memory footprint savings via the use of depth-wise separated convolutions.
• In the object detection space, YoloV3 (Redmon et Farhadi) and Single-Shot Multibox Detector (SSD, Liu et al) are known to be usually faster and less memory-intensive than FasterRCNN (Ren et al), yet also generally less accurate.
• For k-nearest neighbor (k-NN) search problems, a number of approximate k-NN algorithm have been developed and enable fast approximate search: FAISS (related paper: Billion-Scale Similarity Search with GPU, Johnson et al), HNSW (related paper: Efficient and robust approximate nearest neighbor search using Hierarchical Navigable Small World graphs, Malkov et Yashunin), or Annoy

Be aware that unlike hardware and software optimization axis, changing model algorithm can significantly alter accuracy.

4. System Architecture Considerations
Beyond improving the model inference functionality, other architectural concepts can contribute to reducing latency and costs of an end-to-end real-time machine learning prediction pipeline:

• Caching: if you expect a high level of prediction requests with the same input, it is worth considering creating a cache layer in front of the model. A cache can store model input-prediction pairs during an appropriate time-to-live, so that a prediction request over an input that was seen already can be fetched from the cache at a lower latency and cost than re-running the model prediction. Typical relevant options for this layer are fast key-value stores such as Amazon DynamoDB and Redis.
• Locality: localizing the model as close as possible to the consumption point of the prediction may contribute to significant latency improvement. For example, running one model.predict() prediction call on a vector of 127 records with the Scikit-Learn random forest regressor of the public SageMaker Scikit-Learn integration demo notebook takes about 9ms of local execution time on a single-machine ml.c5.large Amazon SageMaker endpoint, but takes 25ms when model serving, serialization/de-serialization and network communication time and are taken into account (from the endpoint to a client instance making the call from within the same region and the same account - both numbers are averages over 10 invocations spaced by multiple seconds). Footprint optimizations mentioned in 1, 2 and 3 above are especially relevant when models are deployed on edge hardware with limited power supply and computing capacity.
• Serverless deployment: For models facing an unpredictable traffic and not requiring GPU a serverless compute engine - for example AWS Lambda - can be a relevant inferencing platform. The pay-per-inference format makes it economical for reasonable amounts of traffic. It also comes with the benefit of simplifying the deployment stack by removing the model server layer. Numerous blog posts show its successful use for inference of deep MXNet models:
• Serving deep learning at Curalate with Apache MXNet, AWS Lambda, and Amazon Elastic Inference
• Build, test, and deploy your Amazon Sagemaker inference models to AWS Lambda

In summary, there are multiple axis to explore to reduce the latency, cost and hardware footprint of your real-time model deployments, and many of those ideas are orthogonal and can be combined. Apache MXNet provides a growing list of tools supporting all flavors of software-based inference improvement, from model architectures to web serving, including compilation and quantization. Please do not hesitate to contribute to those projects and reach out to the community on the forum discuss.mxnet.io/!

Many thanks to Jonathan Taws, Chris Fregly, Will Badr, Tanya Roostaeyan, Madhuri Peri, Frank Liu, Jonathan Chung and Nikolay Ulmasov who provided helpful reviews on the topic

Faster, Cheaper, Leaner: improving real-time ML inference using Apache MXNet was originally published in Apache MXNet on Medium, where people are continuing the conversation by highlighting and responding to this story.

import collections as col

# text stores a piece of text
words = col.Counter(text.split(sep=’ ‘))

Store result in a dictionary

result = dict(col.Counter(text.split(sep=’ ‘))

Store in a pandas data frame

result = pd.DataFrame(col.Counter(text.split(sep=’ ‘)),
                      columns=['word','count'])

The most_common method allows to restrict the counter to top N occurrences

topten = col.Counter(text.split(sep=’ ‘)).most_common(10)

percentiles (p in [0;100])

np.nanpercentile(data[measure], p)

log10

np.log10(data[measure])

power of 10

np.power(10, data[measure])

Finding columns with infinite values (-inf or inf)

def list_inf(df):
    return [f for f in df.columns if np.isinf(df[f]).any()]

Creating a binary from threshold T of metric m

np.sign(data[m] + T)

Weighted average

np.average(df[measure], weights=df[weights])

def compound_growth(L):
    return spst.gmean(1 + np.array(L)) — 1

import operator as op
import pandas as pd
import functools as fct  # reduce is in functools in python 3+

def listagg(L):
    return fct.reduce(op.concat, L)

df.groupby(dim)[textcolumn].apply(listagg)

import networkx as nx

def plot_graph(nodes, edges, labels=False, node_size=False,
               node_color='r', arrows=False, alpha=0.1):

    """nodes: list of nodes. Used as labels if labels=True
       edges: list of edges. Edges are tuples.
       node_size: size or list of sizes
       node_color: color or list of colors
    """

    G = nx.DiGraph()
    G.add_nodes_from(nodes)
    G.add_edges_from(edges)
    nx.draw(G, with_labels=labels, node_color=node_color,
            node_size=node_size, arrows=arrows, alpha=alpha)

Example with a simple fictional dataset:

import pandas as pd

# Node data
Stations = pd.DataFrame([['Station1', 2500],['Station2',1210],
                         ['Station3', 90],['Station4', 312],
                         ['Station5', 11],['Station6', 190]],
                        columns=['Station', 'Traffic'])

# Edge data
routes = pd.DataFrame([['Station1', 'Station2'],
                       ['Station1', 'Station4'],
                       ['Station1', 'Station5'],
                       ['Station2', 'Station3'],
                       ['Station2', 'Station5'],
                       ['Station2', 'Station6'],
                       ['Station3', 'Station5'],
                       ['Station4', 'Station1'],
                       ['Station4', 'Station2'],
                       ['Station4', 'Station6'],
                       ['Station5', 'Station6']],
                      columns=['start', 'stop'])

plot_graph(nodes=Stations['Station'],
           edges=[tuple(row) for row in routes.values],
           labels=True,
           node_color='orange',
           node_size=Stations['Traffic'],
           alpha=0.8)

Example with an anonymized Facebook dataset (Stanford Large Network Dataset Collection); 4,039 nodes 88,233 edges:

data = pd.read_csv('facebook_combined.txt', sep='\t')

# At this stage edges are strings. We transform them to tuples:
edges = [tuple(x[0].split(' ')) for x in data.values]

# We use a set over a nested list comprehension to make a deduplicated list of nodes
nodes = set([item for edge in edges for item in edge])

plot_graph(nodes, edges, alpha=0.05)

company = collections.namedtuple('company', 'ceo, ticker, price')

Create instances from lists

# load some data from a list, e.g.:
data = [‘Jeff Bezos’, ‘AMZN’, 997]

# _make method
Amazon = company._make(data)

# or list unpacking; same result
Amazon = company(*data)

Get data back

Amazon.ticker

Statsmodels:

import statsmodels.api as sm

# fit
X = sm.add_constant(X)
model = sm.OLS(Y,X).fit()

# get an exhaustive model summary
model.summary()

# coefficients
model.params

# rsquared
model.rsquared

# predict out-of-sample
model.predict(Xpred)

Sklearn:

import sklearn

# fit
model = sklearn.linear_model.LinearRegression()
model.fit(X, Y)

# coefficients
model.coef_

# intercept
model.intercept_

# sklearn does not have rsquare as a model attribute

# predict out-of-sample
pred = model.predict(Xpred)

Numpy:

import numpy as np

np.linalg.lstsq(x, y)

Scipy:

import scipy as sp

# returns parameters in one method
slope, intercept, r_value, p_value, std_err = \
    sp.stats.linregress(x,y)

Values in list:

df[df.col.isin(L)]

String like:

df[df.col.str.contains('text')]

Not null:

df[~df.col.isnull()]  # using filter negation ~