Blaze: Next Generation Numpy

Blaze is a next-generation NumPy sponsored by Continuum Analytics. It is designed as a foundational set of abstractions on which to build out-of-core and distributed algorithms. Blaze generalizes many of the ideas found in popular PyData projects such as Numpy, Pandas, and Theano into one generalized data-structure. Together with a powerful array-oriented virtual machine and run-time, Blaze will be capable of performing efficient linear algebra and indexing operations on top of a wide variety of data backends.

1. . . . . . . . . . . . . . Blaze Stephen Diehl Continuum Analytics Stephen Diehl Blaze
2. . . . . . . . . Warning Extreme talk, code and math ahead! Stephen Diehl Blaze
3. . . . . . . . . Problems with NumPy • Limited support for heterogeneous data • Missing values ( na ) • Uses a single buffer with strides overlayed • Cannot resize in place • Labeled arrays • No support for arrays larger than available memory ( out-of-core ) • No support for Distributed Arrays Stephen Diehl Blaze
4. . . . . . . . . What is Blaze How do we extend NumPy so that we can map array computations onto modern hardware? • Blaze is • a generalization of NumPy • a datashape description language • a generalized notion of data access • a way to view data stores as arrays and tables Stephen Diehl Blaze
5. . . . . . . . . Inspirations Object Oriented Programming Object 1 Object 3 Object 2 Attr 1 Attr 2 Attr 3 Attr 1 Attr 2 Attr 3 Attr 1 Attr 2 Attr 3 Stephen Diehl Blaze
6. . . . . . . Array Oriented Programming Object 1 Object 2 Object 3 Attr 1 Attr 2 Attr 3 Stephen Diehl Blaze
7. . . . . . . . . Array Oriented Programming • APL • First class arrays • Chapel / ZPL • Domain Maps • Tensor Expressions Stephen Diehl Blaze
8. . . . . . . . . APL Dialects e notion of array shape is extremely fluid. > 3 4 2 $ a 0 1 2 3 4 5 6 7 8 0 1 2 3 4 5 6 7 8 0 1 2 3 4 5 ... Stephen Diehl Blaze
9. . . . . . . forall ijkl in A do A[ijkl] = reduce f [ijl] + reduce g [jkl] Very low barrier between the high-level mathematical notation and the programming notation. (a × b)i = εijkajbk Stephen Diehl Blaze
10. . . . . . . . Functional Programming . . . . . . . . • Haskell • Types • Lazy Evaluation • Data Parallel Haskell • Stream Fusion • Repa Stephen Diehl Blaze
11. . . . . . . . . End-User Perspectives Blaze is built around separation of concerns. • Domain experts • Richer structure to express high level ideas. • Algorithm writers • More information ( type, shape, layout ) to do clever optimizations. • Researchers • A platform in which to explore data and task parallelism. . Core Idea . . . . . . . . Leverage the existing community and projects. Stephen Diehl Blaze
12. . . . . . . . . Ecosystem Requirements Numpy DataFrame Business Data • Matrices • Multidimensional Models • Data Cubes • OLAP / MDX • Timeseries • Hierarchical Scalar Categorical Dimensions Structures Vectors Sets Sequential Random Data Access Strong Coupling of Values and Data Locality Weak Coupling of Values and Data Locality Walk kernels over grid structure of the array Statistical Aggregations Stephen Diehl Blaze
13. . . . . . . . . Current Work • Data Structures • PyTables • carray • varray • ctable • pandas • memmap, masked array • Computation • theano • numexpr • clyther Stephen Diehl Blaze
14. . . . . . . . . NumPy • What defines a NumPy array? In[0]: A = np.array([[1, 2], [3, 4]], dtype=np.int32) array([[1, 2], [3, 4]], dtype=int32) Stephen Diehl Blaze
15. . . . . . . . . Data In[1]: bytes(A.data) Out[1]: b'\x01\x00\x00\x00\x02\x00\x00\x00\x03\x00\x00...' Stephen Diehl Blaze
16. . . . . . . . . Shape In[2]: A.shape Out[2]: (2,2) Stephen Diehl Blaze
17. . . . . . . . . DType In[3]: A.dtype Out[3]: dtype('int32') Stephen Diehl Blaze
18. . . . . . . . . Strides In[4]: A.strides Out[4]: (8,4) Stephen Diehl Blaze
19. . . . . . . . . Flags In[5]: A.flags Out[5]: C_CONTIGUOUS : True F_CONTIGUOUS : False OWNDATA : True WRITEABLE : True ALIGNED : True UPDATEIFCOPY : False Stephen Diehl Blaze
20. . . . . . . How do we make extend NumPy so that we can map array computations onto modern hardware? • How do we generalize dtype system? • How do we construct Numpy arrays out of multiple buffers? • How do we dispatch Numpy expressions to specialized kernels? Stephen Diehl Blaze
21. . . . . . . . . Datashapes How do we make dtype more general? Stephen Diehl Blaze
22. . . . . . . . . Datashape DType Shape Datashape Stephen Diehl Blaze
23. . . . . . . . Datashape . . . . . . . . Datashapes are an extension of dtype, expressed as a small DSL. Stephen Diehl Blaze
24. . . . . . . . . Directory of Images Endow existing data with structure. $ ls /folder1 lena1.png lena2.png /folder2 lena3.png lena4.png /folder3 lena5.png lena6.png ... /folderN Stephen Diehl Blaze
25. . . . . . . N M 800 R G B A 0 600 (N, M, 800, 600, RGBA) Stephen Diehl Blaze
26. . . . . . . . . CSVs Endow existing data with structure. $ ls /particle_accelerator measurement1.csv measurement2.csv measurement3.csv Stephen Diehl Blaze
27. . . . . . . • Multiple views of the same hypercube. • Fluid interpretation of shape. Ax Ay Az Bx By Bz 0 3.1 4.1 5.9 2.6 5.3 5.8 1 2.7 1.8 2.8 1.8 2.8 4.5 2 0.5 7.7 2.1 5.6 6.4 9.0 Stephen Diehl Blaze
28. . . . . . . (A 0 , B 0 ) (A 2 , B 2 ) (A 1 , B 1 ) (A 3 , B 3 ) (A 5 , B 5 ) (A 4 , B 4 ) (A 6 , B 6 ) (A 8 , B 8 ) (A 7 , B 7 ) y x A y A x B x B y 3 1 0 1 2 7 1 5 4 9 6 4 7 8 3 1 1 2 0 8 1 5 8 2 Stephen Diehl Blaze
29. . . . . . . 5, 5, int32 N, M, 800, 600, RGBA 300, Record(uid=int32, name=string) Var(25), Record(uid=int32, name=string) 200, 300, Either(int32, na) 10, 10, Quaternion Stephen Diehl Blaze
30. . . . . . . . . Record class Stock(RecordClass): name = string open = decimal close = decimal max = int64 min = int64 @derived def mid(self): # Types are inferred return (self.min + self.max)/2 Stephen Diehl Blaze
31. . . . . . . . . Type Constructors • Alias types RGBA = Record(r=int32, g=int32, b=int32, a=int8) Stock = Record(open=decimal, close=decimal, name=string) • Parametric types SquareMatrix A = A, A Mesh A B = 5, 5, Record(x=A, y=B) Portfolio R = R, Stock Stephen Diehl Blaze
32. . . . . . . . . Usage In[0]: stream = blaze.open('yahoo://GOOG:20120100') In[1]: nd = NDTable(stream, datashape='Var(100), Stock') In[2]: nd[0] name open close max min mid Google Inc. 681.79 705.92 695 621 658 Syntax not final. Stephen Diehl Blaze
33. . . . . . . . . Coordinates • How do we construct NumPy arrays out of multiple buffers? Stephen Diehl Blaze
34. . . . . . . . . Memory Access void * PyArray_GetPtr(PyArrayObject *obj, npy_intp* ind) { int n = obj->nd; npy_intp *strides = obj->strides; char *dptr = obj->data; while (n--) { // Core pointer arithmetic of Numpy dptr += (*strides++) * (*ind++); } return (void *)dptr; } Stephen Diehl Blaze
35. . . . . . . . . Numpy Data Access Sn = itemsize × m ∏ n+1 shapen I = (i, j, k, ...) Ai,j,k = S · I Stephen Diehl Blaze
36. . . . . . . . . Blaze Data Access Sn = m ∏ n+1 datashapen I = (i, j, k, ...) Ai,j,k = f(S; I) Stephen Diehl Blaze
37. . . . . . . • Stride vectors are now functions of indexers and the datashape. Sn = m ∏ n+1 datashapen • Generalized indexers: scalars, labels, references. I = (i, j, k, ...) Ai,j,k = f(S; I) • In the case of strided memory access: f = (·) . Stephen Diehl Blaze
38. . . . . . . . . Coordinates • Coordinates • Atlas of regions in memory to resolve partitions. • Charts of memory with coordinate systems. • Coordinate transformations between charts. • Index space transformations are homomorphisms between charts. Stephen Diehl Blaze
39. . . . . . . 0 0 0 0 0 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 1 1 1 1 vstack ( ( , 3 4 x 3 4 x T 0 0 0 0 0 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 1 1 1 1 6 4 x 4 6 x 0 0 0 0 0 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 1 1 1 1 f = vstack p = (0,3) g = transpose Stephen Diehl Blaze
40. . . . . . . f(i, j) = (i, j) + p g(i, j) = (j, i) (f.g)(i, j) = (j, i) + g(p) Stephen Diehl Blaze
41. . . . . . . (3,5) 4 6 x 0 0 0 0 0 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 1 1 1 1 Stephen Diehl Blaze
42. . . . . . . (2,3) 0 0 0 0 0 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 1 1 1 1 Stephen Diehl Blaze
43. . . . . . . (3,2) 0 0 0 0 0 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 1 1 1 1 Stephen Diehl Blaze
44. . . . . . . \x00\x00\x00\x00\x00\x00\xf0?\x00\x00\x00\x00\x00\x00 \xf0?\x00\x00\x00\x00\x00\x00\xf0?\x00\x00\x00\x00\x00 \x00\xf0?\x00\x00\x00\x00\x00\x00\xf0?\x00\x00\x00\x00 \x00\x00\xf0?\x00\x00\x00\x00\x00\x00\xf0?\x00\x00\x00 \x00\x00\x00\xf0?\x00\x00\x00\x00\x00\x00\xf0?\x00\x00 \x00\x00\x00\x00\xf0?\x00\x00\x00\x00\x00\x00\xf0?\x00 \x00\x00\x00\x00\x00\xf0?' int32 = 4 bytes offset = 16 0 0 0 0 0 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 1 1 1 1 Stephen Diehl Blaze
45. . . . . . . . . Blaze Data Descriptors At some level our coordinates manifest simply as byte offsets in some storage medium. • Memory Like • Arbitrary Slices • Random Seeks • File Like • Chunks • Random Seeks • Stream Like • Chunks • Sequential Seeks Stephen Diehl Blaze
46. . . . . . . Chunked Streaming Contigious Copy Seek Stephen Diehl Blaze
47. . . . . . . . . Dispatch How do we dispatch NumPy expressions to specialized kernels? Stephen Diehl Blaze
48. . . . . . . . . Indexables . Lazy Evaluation . . . . . . . . Expressions are not evaluated until forced. • Deferred table. • Appends graph nodes onto existing graph. • Data processing scripts, fine grained control over evaluation. • Immediate table. • Forces graph immediately. • Data exploration at REPL, when you just want a number. Stephen Diehl Blaze
49. . . . . . . NDTable NDArray Deferred Indexable Immediate Deferred Immediate Stephen Diehl Blaze
50. . . . . . . . . Expression Graph In[0]: x = a + b In[1]: type(x) <type 'blaze.graph.ArrayNode'> Stephen Diehl Blaze
51. . . . . . . . . Expression Graph In[0]: x = a + b In[1]: type(x) <type 'blaze.graph.ArrayNode'> In[2]: e = LA.eig(x) In[3]: e.eval() Out[3]: array([[ 0.70710678+0.j , 0.70710678+0.j ], [ 0.00000000-0.70710678j, 0.00000000+0.70710678j]]) Stephen Diehl Blaze
52. . . . . . . . . Expression Typing Typing and polymorphic operators coexist. class add(BinaryOp): signature = 'a -> b -> b' dom = [scalar, array] In[0]: expr = 2 + [1,2,3] In[1]: expr = [1,2,3] + 2 In[2]: expr.dom [int32[:], int32[:]] In[3]: expr.cod int32[:] Stephen Diehl Blaze
53. . . . . . . . . Transformations f(a, b) :: (Indexable, Indexable) → Indexable • Index space transformation • How coordinates transform under f. • Value space transformation • How values transform under f. • Metadata space transformation • How metadata is merged from operands a, b. Preserving type and metadata information allows us to write more clever algorithms that avoid needless seeks and scans. Stephen Diehl Blaze
54. . . . . . . . . Nuances of Evaluation • Forcing Evaluation • Many properties of arrays are not knowable a priori. Evaluation is forced when we need to peek into a opaque type. • Can't write coordinate transformations in terms of coordinates we can't describe. • Side effectful operations. ( numpy.fill, numpy.put ) • (reduce, groupby) also force evaluation. In general operations which are not index space homomorphisms. • Just like NumPy. ere are benefits to writing in functional style, but one can still drop into raw memory manipulation when needed. Stephen Diehl Blaze
55. . . . . . . . . Type Unification Dtype Promotion Broadcasting Type Unification Stephen Diehl Blaze
56. . . . . . . NumPy A 5 x 4 B 1 broadcast(A,B) 5 x 4 Blaze A 3, 2, Var(10), int32 B 1, Var(100), float32 unify(A,B) 3, 2, Var(100), float32 Stephen Diehl Blaze
57. . . . . . . . . Metadata • Data Warehousing • Provenance: Where does this number come from? • Propagate metadata about source through ETL pipeline. • Orthogonal matrix • A−1 → AT • Value Space Transform → Index Space Transform • Matrix density • density(A) + density(B) in terms of density(A,B) Provide an extensible algebra of propagating information through the expression graph. Stephen Diehl Blaze
58. . . . . . . . . Global Arrays Blaze as the lingua franca of Python data exchange. • Remote Computation • Push code to data • Remote Data Retrieval • Pull data to code Stephen Diehl Blaze
59. . . . . . . . . Compiler Backend • Code generation • Numba & Minivect ( target LLVM ) • Array VM ( target custom bytecode ) Kernel Specialization Metadata Data Locality Data Layout OpenMP BLAS/LAPACK Cuda Stephen Diehl Blaze
60. . . . . . . . . Blaze Blaze is still a work in progress. Look for more details at the end of November. • Blaze will be open source! • Blaze will have tight integration with Numba. • Blaze Studio will have premium features on top of Blaze and Numba and will include IOPro. Stephen Diehl Blaze
61. . . . . . . . . Blaze Project . People . . . . . . . . • Travis Oliphant • Peter Wang • Mark Wiebe ( numpy ) • Francesc Alted ( pytables, numexpr ) • Mark Florisson ( minivect, numba ) • Stephen Diehl . Contact . . . . . . . . • Github: github.com/ContinuumIO/blaze • URL: http://www.continuum.io • Twitter: @ContinuumIO Stephen Diehl Blaze