Array programming with NumPy — TEI-tagged view

Charles R. HarrispersName, K. Jarrod MillmanpersName, Stéfan J. van der WaltpersName, Ralf GommerspersName, Pauli VirtanenpersName, David CournapeaupersName, Eric WieserpersName, Julian TaylorpersName, Sebastian BergpersName, Nathaniel J. SmithpersName, Robert KernpersName, Matti PicuspersName, Stephan HoyerpersName, Marten H. van KerkwijkpersName, Matthew BrettpersName, Allan HaldanepersName, Jaime Fernández del RíopersName, Mark WiebepersName, Pearu PetersonpersName, Pierre Gérard-MarchantpersName, Kevin SheppardpersName, Tyler ReddypersName, Warren WeckesserpersName, Hameer AbbasipersName, Christoph GohlkepersName, Travis E. OliphantpersName (2020date). “Array programming with NumPy.” Naturetitle·j.

Object of study — NumPyterm's array concepts, scientific-Python role, interoperability protocolsterm, development model, and future constraints

Research question — how fundamental array concepts and a shared APIterm organize scientific computation while allowing specialized arrayterm implementations to interoperate

Contribution — a review of NumPyterm's array-programming model, ecosystem position, and transition into the standard coordinating APIterm for Pythonterm array technologies

Full text — tagged

Array programmingterm with NumPyterm Abstract Array programmingterm provides a powerful, compact and expressive syntax for accessing, manipulating and operating on data in vectors, matrices and higher-dimensional arrays. NumPyterm is the primary array programmingterm library for the Pythonterm language. It has an essential role in research analysis pipelines in fields as diverse as physics, chemistry, astronomy, geoscience, biology, psychology, materials science, engineering, finance and economics. For example, in astronomy, NumPyterm was an important part of the software stack used in the discovery of gravitational waves 1 and in the first imaging of a black hole 2 . Here we review how a few fundamental array concepts lead to a simple and powerful programming paradigm for organizing, exploring and analysing scientific data. NumPyterm is the foundation upon which the scientific Python ecosystemterm is constructed. It is so pervasive that several projects, targeting audiences with specialized needs, have developed their own NumPy-like interfaces and array objects. Owing to its central position in the ecosystem, NumPyterm increasingly acts as an interoperability layerterm between such array computation libraries and, together with its application programming interfaceterm (APIterm), provides a flexible framework to support the next decade of scientific and industrial analysis.

Subject terms: Computational neuroscience, Solar physics, Computational science, Computer science, Software Main Two Pythonterm array packages existed before NumPyterm. The Numericterm package was developed in the mid-1990s and provided array objects and array-aware functionsterm in Pythonterm. It was written in C and linked to standard fast implementations of linear algebra 3 , 4 . One of its earliest uses was to steer C++ applications for inertial confinement fusion research at Lawrence Livermore National Laboratory 5 . To handle large astronomical images coming from the Hubble Space Telescope, a reimplementation of Numericterm, called Numarrayterm, added support for structured arrays, flexible indexingterm, memory mapping, byte-order variants, more efficient memory use, flexible IEEE 754-standard error-handling capabilities, and better type-casting rules 6 . Although Numarrayterm was highly compatible with Numericterm, the two packages had enough differences that it divided the community; however, in 2005 NumPyterm emerged as a ‘best of both worldsterm’ unification 7 —combining the features of Numarrayterm with the small-array performance of Numericterm and its rich C APIterm.

Now, 15 years later, NumPyterm underpins almost every Pythonterm library that does scientific or numerical computation 8 – 11 , including SciPyterm 12 , Matplotlibterm 13 , pandas 14 , scikit-learn 15 and scikit-image 16 . NumPyterm is a community-developed, open-source library, which provides a multidimensional Pythonterm array object along with array-aware functionsterm that operate on it. Because of its inherent simplicity, the NumPy arrayterm is the de facto exchange format for array data in Pythonterm.

NumPyterm operates on in-memory arrays using the central processing unit (CPU). To utilize modern, specialized storage and hardware, there has been a recent proliferation of Pythonterm array packages. Unlike with the Numarrayterm–Numericterm divide, it is now much harder for these new libraries to fracture the user community—given how much work is already built on top of NumPyterm. However, to provide the community with access to new and exploratory technologies, NumPyterm is transitioning into a central coordinating mechanismterm that specifies a well defined array programmingterm APIterm and dispatches it, as appropriate, to specialized arrayterm implementations.

NumPy arraysterm The NumPy arrayterm is a data structureterm that efficiently stores and accesses multidimensional arraysterm 17 (also known as tensorsterm), and enables a wide variety of scientific computation. It consists of a pointer to memoryterm, along with metadataterm used to interpret the data stored there, notably ‘data typeterm’, ‘shapeterm’ and ‘stridesterm’ (Fig. 1a ).

The NumPy arrayterm incorporates several fundamental array concepts. a , The NumPy arrayterm data structureterm and its associated metadataterm fields. b , Indexingterm an array with slices and steps. These operations return a ‘viewterm’ of the original data. c , Indexingterm an array with masks, scalar coordinates or other arrays, so that it returns a ‘copy’ of the original data. In the bottom example, an array is indexed with other arrays; this broadcasts the indexingterm arguments before performing the lookup. d , Vectorizationterm efficiently applies operations to groups of elements. e , Broadcastingterm in the multiplication of two-dimensional arrays. f , Reductionterm operations act along one or more axesterm·var. In this example, an array is summed along select axesterm·var to produce a vector, or along two axesterm·var consecutively to produce a scalar. g , Example NumPyterm code, illustrating some of these concepts.

The data typeterm describes the nature of elements stored in an array. An array has a single data typeterm, and each element of an array occupies the same number of bytes in memory. Examples of data typesterm include real and complex numbers (of lower and higher precision), strings, timestamps and pointers to Pythonterm objects.

The shapeterm of an array determines the number of elements along each axis, and the number of axesterm·var is the dimensionalityterm of the array. For example, a vector of numbers can be stored as a one-dimensional array of shapeterm N , whereas colour videos are four-dimensional arrays of shapeterm ( T , M , N , 3).

Stridesterm are necessary to interpret computer memory, which stores elements linearly, as multidimensional arraysterm. They describe the number of bytes to move forward in memory to jump from row to row, column to column, and so forth. Consider, for example, a two-dimensional array of floating-point numbers with shapeterm (4, 3), where each element occupies 8 bytes in memory. To move between consecutive columns, we need to jump forward 8 bytes in memory, and to access the next row, 3 × 8 = 24 bytes. The stridesterm of that array are therefore (24, 8). NumPyterm can store arrays in either C or Fortran memory orderterm, iterating first over either rows or columns. This allows external libraries written in those languages to access NumPy arrayterm data in memory directly.

Users interact with NumPy arraysterm using ‘indexingterm’ (to access subarrays or individual elements), ‘operatorsterm’ (for example, +, − and × for vectorized operations and @ for matrix multiplication), as well as ‘array-aware functionsterm’; together, these provide an easily readable, expressive, high-level APIterm for array programmingterm while NumPyterm deals with the underlying mechanics of making operations fast.

Indexingterm an array returns single elements, subarrays or elements that satisfy a specific condition (Fig. 1b ). Arrays can even be indexed using other arrays (Fig. 1c ). Wherever possible, indexingterm that retrieves a subarray returns a ‘viewterm’ on the original array such that data are shared between the two arrays. This provides a powerful way to operate on subsets of array data while limiting memory usage.

To complement the array syntax, NumPyterm includes functions that perform vectorized calculations on arrays, including arithmetic, statistics and trigonometry (Fig. 1d ). Vectorizationterm—operating on entire arrays rather than their individual elements—is essential to array programmingterm. This means that operations that would take many tens of lines to express in languages such as C can often be implemented as a single, clear Pythonterm expression. This results in concise code and frees users to focus on the details of their analysis, while NumPyterm handles looping over array elements near-optimally—for example, taking stridesterm into consideration to best utilize the computer’s fast cache memory.

When performing a vectorized operation (such as addition) on two arrays with the same shapeterm, it is clear what should happen. Through ‘broadcastingterm’ NumPyterm allows the dimensions to differ, and produces results that appeal to intuition. A trivial example is the addition of a scalar value to an array, but broadcastingterm also generalizes to more complex examples such as scaling each column of an array or generating a grid of coordinates. In broadcastingterm, one or both arrays are virtually duplicated (that is, without copying any data in memory), so that the shapesterm of the operands match (Fig. 1d ). Broadcastingterm is also applied when an array is indexed using arrays of indices (Fig. 1c ).

Other array-aware functionsterm, such as sum, mean and maximum, perform element-by-element ‘reductionsterm’, aggregating results across one, multiple or all axesterm·var of a single array. For example, summing an n -dimensional array over d axesterm·var results in an array of dimension n − d (Fig. 1f ).

NumPyterm also includes array-aware functionsterm for creating, reshaping, concatenating and padding arrays; searching, sorting and counting data; and reading and writing files. It provides extensive support for generating pseudorandom numbersterm, includes an assortment of probability distributions, and performs accelerated linear algebra, using one of several backends such as OpenBLASterm 18 , 19 or Intel MKLterm optimized for the CPUs at hand (see Supplementary Methods for more details).

Altogether, the combination of a simple in-memory array representation, a syntax that closely mimics mathematics, and a variety of array-aware utility functions forms a productive and powerfully expressive array programmingterm language.

Scientific Python ecosystemterm Pythonterm is an open-source, general-purpose interpreted programming language well suited to standard programming tasks such as cleaning data, interacting with web resources and parsing text. Adding fast array operations and linear algebra enables scientists to do all their work within a single programming language—one that has the advantage of being famously easy to learn and teach, as witnessed by its adoption as a primary learning language in many universities.

Even though NumPyterm is not part of Pythonterm’s standard library, it benefits from a good relationship with the Pythonterm developers. Over the years, the Pythonterm language has added new features and special syntax so that NumPyterm would have a more succinct and easier-to-read array notation. However, because it is not part of the standard library, NumPyterm is able to dictate its own release policies and development patterns.

SciPyterm and Matplotlibterm are tightly coupled with NumPyterm in terms of history, development and use. SciPyterm provides fundamental algorithms for scientific computing, including mathematical, scientific and engineering routines. Matplotlibterm generates publication-ready figures and visualizations. The combination of NumPyterm, SciPyterm and Matplotlibterm, together with an advanced interactive environment such as IPythonterm 20 or Jupyterterm 21 , provides a solid foundation for array programmingterm in Pythonterm. The scientific Python ecosystemterm (Fig. 2 ) builds on top of this foundation to provide several, widely used technique-specific libraries 15 , 16 , 22 , that in turn underlie numerous domain-specific projects 23 – 28 . NumPyterm, at the base of the ecosystem of array-aware libraries, sets documentation standards, provides array testing infrastructure and adds build support for Fortran and other compilers.

NumPyterm is the base of the scientific Python ecosystemterm. Essential libraries and projects that depend on NumPyterm’s APIterm gain access to new array implementations that support NumPyterm’s array protocolsterm·def (Fig. 3 ).

Many research groups have designed large, complex scientific libraries that add application-specific functionality to the ecosystem. For example, the eht-imaging library 29 , developed by the Event Horizon Telescope collaboration for radio interferometry imaging, analysis and simulation, relies on many lower-level components of the scientific Python ecosystemterm. In particular, the EHT collaboration used this library for the first imaging of a black hole. Within eht-imaging, NumPy arraysterm are used to store and manipulate numerical data at every step in the processing chain: from raw data through calibration and image reconstruction. SciPyterm supplies tools for general image-processing tasks such as filtering and image alignment, and scikit-image, an image-processing library that extends SciPyterm, provides higher-level functionality such as edge filters and Hough transforms. The ‘scipyterm.optimize’ module performs mathematical optimization. NetworkX 22 , a package for complex network analysis, is used to verify image comparison consistency. Astropy 23 , 24 handles standard astronomical file formats and computes time–coordinate transformations. Matplotlibterm is used to visualize data and to generate the final image of the black hole.

The interactive environment created by the array programmingterm foundation and the surrounding ecosystem of tools—inside of IPythonterm or Jupyterterm—is ideally suited to exploratory data analysisterm. Users can fluidly inspect, manipulate and visualize their data, and rapidly iterate to refine programming statements. These statements are then stitched together into imperative or functional programs, or notebooks containing both computation and narrative. Scientific computing beyond exploratory work is often done in a text editor or an integrated development environment (IDE) such as Spyder. This rich and productive environment has made Pythonterm popular for scientific research.

To complement this facility for exploratory work and rapid prototyping, NumPyterm has developed a culture of using time-tested software engineering practices to improve collaboration and reduce error 30 . This culture is not only adopted by leaders in the project but also enthusiastically taught to newcomers. The NumPyterm team was early to adopt distributed revision controlterm and code reviewterm to improve collaboration on code, and continuous testingterm that runs an extensive battery of automated tests for every proposed change to NumPyterm. The project also has comprehensive, high-quality documentation, integrated with the source code 31 – 33 .

This culture of using best practices for producing reliable scientific software has been adopted by the ecosystem of libraries that build on NumPyterm. For example, in a recent award given by the Royal Astronomical Society to Astropy, they state: “The Astropy Project has provided hundreds of junior scientists with experience in professional-standard software development practices including use of version control, unit testing, code reviewterm and issue tracking procedures. This is a vital skill set for modern researchers that is often missing from formal university education in physics or astronomy” 34 . Community members explicitly work to address this lack of formal education through courses and workshops 35 – 37 .

The recent rapid growth of data science, machine learning and artificial intelligence has further and dramatically boosted the scientific use of Pythonterm. Examples of its important applications, such as the eht-imaging library, now exist in almost every discipline in the natural and social sciences. These tools have become the primary software environment in many fields. NumPyterm and its ecosystem are commonly taught in university courses, boot camps and summer schools, and are the focus of community conferences and workshops worldwide. NumPyterm and its APIterm have become truly ubiquitous.

Array proliferation and interoperability NumPyterm provides in-memory, multidimensional, homogeneously typed (that is, single-pointer and strided) arrays on CPUs. It runs on machines ranging from embedded devices to the world’s largest supercomputers, with performance approaching that of compiled languages. For most its existence, NumPyterm addressed the vast majority of array computation use cases.

However, scientific datasets now routinely exceed the memory capacity of a single machine and may be stored on multiple machines or in the cloud. In addition, the recent need to accelerate deep-learning and artificial intelligence applications has led to the emergence of specialized accelerator hardware, including graphics processing unitsterm (GPUs), tensor processing unitsterm (TPUs) and field-programmable gate arraysterm (FPGAs). Owing to its in-memory data modelterm, NumPyterm is currently unable to directly utilize such storage and specialized hardwareterm·var. However, both distributed data and also the parallel execution of GPUs, TPUs and FPGAs map well to the paradigm of array programmingterm: therefore leading to a gap between available modern hardware architectures and the tools necessary to leverage their computational power.

The community’s efforts to fill this gap led to a proliferation of new array implementations. For example, each deep-learning framework created its own arrays; the PyTorchterm 38 , Tensorflowterm 39 , Apache MXNet 40 and JAXterm arrays all have the capability to run on CPUs and GPUs in a distributed fashion, using lazy evaluation to allow for additional performance optimizations. SciPyterm and PyData/Sparse both provide sparse arraysterm, which typically contain few non-zero valuesterm and store only those in memory for efficiency. In addition, there are projects that build on NumPy arraysterm as data containers, and extend its capabilities. Distributed arraysterm are made possible that way by Daskterm, and labelled arraysterm—referring to dimensions of an array by name rather than by index for clarity, compare x[:, 1] versus x.loc[:, 'time']—by xarrayterm 41 .

Such libraries often mimic the NumPyterm APIterm, because this lowers the barrier to entry for newcomers and provides the wider community with a stable array programming interfaceterm. This, in turn, prevents disruptive schisms such as the divergence between Numericterm and Numarrayterm. But exploring new ways of working with arrays is experimental by nature and, in fact, several promising libraries (such as Theano and Caffe) have already ceased development. And each time that a user decides to try a new technology, they must change import statements and ensure that the new library implements all the parts of the NumPyterm APIterm they currently use.

Ideally, operating on specialized arraysterm using NumPyterm functions or semantics would simply work, so that users could write code once, and would then benefit from switching between NumPy arraysterm, GPU arrays, distributed arraysterm and so forth as appropriate. To support array operations between external array objects, NumPyterm therefore added the capability to act as a central coordination mechanism with a well specified APIterm (Fig. 2 ).

To facilitate this interoperability, NumPyterm provides ‘protocolsterm’ (or contracts of operationterm), that allow for specialized arraysterm to be passed to NumPyterm functions (Fig. 3 ). NumPyterm, in turn, dispatches operationsterm to the originating library, as required. Over four hundred of the most popular NumPyterm functions are supported. The protocolsterm are implemented by widely used libraries such as Daskterm, CuPy, xarrayterm and PyData/Sparse. Thanks to these developments, users can now, for example, scale their computation from a single machine to distributed systems using Daskterm. The protocolsterm also compose well, allowing users to redeploy NumPyterm code at scale on distributed, multi-GPU systems via, for instance, CuPy arraysterm embedded in Daskterm arrays. Using NumPyterm’s high-level APIterm, users can leverage highly parallel code execution on multiple systems with millions of cores, all with minimal code changes 42 .

NumPyterm’s APIterm and array protocolsterm·def expose new arrays to the ecosystem. In this example, NumPyterm’s ‘mean’ function is called on a Daskterm array. The call succeeds by dispatching to the appropriate library implementation (in this case, Daskterm) and results in a new Daskterm array. Compare this code to the example code in Fig. 1g .

These array protocolsterm·def are now a key feature of NumPyterm, and are expected to only increase in importance. The NumPyterm developers—many of whom are authors of this Review—iteratively refine and add protocolterm designs to improve utility and simplify adoption.

Discussion NumPyterm combines the expressive power of array programmingterm, the performance of C, and the readability, usability and versatility of Pythonterm in a mature, well tested, well documented and community-developed library. Libraries in the scientific Python ecosystemterm provide fast implementations of most important algorithms. Where extreme optimization is warranted, compiled languages can be used, such as Cython 43 , Numba 44 and Pythran 45 ; these languages extend Pythonterm and transparently accelerate bottlenecks. Owing to NumPyterm’s simple memory model, it is easy to write low-level, hand-optimized code, usually in C or Fortran, to manipulate NumPy arraysterm and pass them back to Pythonterm. Furthermore, using array protocolsterm·def, it is possible to utilize the full spectrum of specialized hardwareterm·var acceleration with minimal changes to existing code.

NumPyterm was initially developed by students, faculty and researchers to provide an advanced, open-source array programmingterm library for Pythonterm, which was free to use and unencumbered by license servers and software protection dongles. There was a sense of building something consequential together for the benefit of many others. Participating in such an endeavour, within a welcoming community of like-minded individuals, held a powerful attraction for many early contributors.

These user–developers frequently had to write code from scratch to solve their own or their colleagues’ problems—often in low-level languages that preceded Pythonterm, such as Fortran 46 and C. To them, the advantages of an interactive, high-level array library were evident. The design of this new tool was informed by other powerful interactive programming languages for scientific computing such as Basis 47 – 50 , Yorick 51 , R 52 and APL 53 , as well as commercial languages and environments such as IDL (Interactive Data Language) and MATLAB.

What began as an attempt to add an array object to Pythonterm became the foundation of a vibrant ecosystem of tools. Now, a large amount of scientific work depends on NumPyterm being correct, fast and stable. It is no longer a small community project, but core scientific infrastructureterm.

The developer culture has matured: although initial development was highly informal, NumPyterm now has a roadmap and a process for proposing and discussing large changes. The project has formal governance structuresterm and is fiscally sponsored by NumFOCUS, a nonprofit that promotes open practices in research, data and scientific computing. Over the past few years, the project attracted its first funded development, sponsored by the Moore and Sloan Foundations, and received an award as part of the Chan Zuckerberg Initiative’s Essentials of Open Source Software programme. With this funding, the project was (and is) able to have sustained focus over multiple months to implement substantial new features and improvements. That said, the development of NumPyterm still depends heavily on contributions made by graduate students and researchers in their free time (see Supplementary Methods for more details).

NumPyterm is no longer merely the foundational array library underlying the scientific Python ecosystemterm, but it has become the standard API for tensor computationterm and a central coordinating mechanismterm between array types and technologies in Pythonterm. Work continues to expand on and improve these interoperability features.

Over the next decade, NumPyterm developers will face several challenges. New devices will be developed, and existing specialized hardwareterm·var will evolve to meet diminishing returns on Moore’s law. There will be more, and a wider variety of, data science practitioners, a large proportion of whom will use NumPyterm. The scale of scientific data gathering will continue to increase, with the adoption of devices and instruments such as light-sheet microscopes and the Large Synoptic Survey Telescope (LSST) 54 . New generation languages, interpreters and compilers, such as Rust 55 , Julia 56 and LLVM 57 , will create new concepts and data structuresterm, and determine their viability.

Through the mechanisms described in this Review, NumPyterm is poised to embrace such a changing landscape, and to continue playing a leading part in interactive scientific computation, although to do so will require sustained fundingterm from government, academia and industry. But, importantly, for NumPyterm to meet the needs of the next decade of data science, it will also need a new generation of graduate students and community contributorsterm to drive it forward.

Supplementary information This file contains Supplementary Methods, including Supplementary Figure 1 and additional references.

Acknowledgements We thank R. Barnowski, P. Dubois, M. Eickenberg, and P. Greenfield, who suggested text and provided helpful feedback on the manuscript. K.J.M. and S.J.v.d.W. were funded in part by the Gordon and Betty Moore Foundation through grant GBMF3834 and by the Alfred P. Sloan Foundation through grant 2013-10-27 to the University of California, Berkeley. S.J.v.d.W., S.B., M.P. and W.W. were funded in part by the Gordon and Betty Moore Foundation through grant GBMF5447 and by the Alfred P. Sloan Foundation through grant G-2017-9960 to the University of California, Berkeley.

Author contributions K.J.M. and S.J.v.d.W. composed the manuscript with input from others. S.B., R.G., K.S., W.W., M.B. and T.R. contributed text. All authors contributed substantial code, documentation and/or expertise to the NumPyterm project. All authors reviewed the manuscript.

Competing interests The authors declare no competing interests.

Footnotes Contributor Information K. Jarrod MillmanpersName, Email: millmanpersName@berkeley.edu.

Stéfan J. van der WaltpersName, Email: stefanv@berkeley.edu Ralf GommerspersName, Email: ralf.gommerspersName@gmail.com.

Supplementary information is available for this paper at 10.1038/s41586-020-2649-2.

References Associated Data Supplementary Materials This file contains Supplementary Methods, including Supplementary Figure 1 and additional references.