NumPy のブロードキャスティングを PHP で模倣する
NumPy の配列には、ブロードキャスティングという仕組みがあります。これは、大きさの異なる多次元配列同士の算術演算に関する仕組みで、要素数が 1 の次元があった場合にそれを拡張し、二つの配列の要素数を揃えて演算結果を得るというものです。NumPy のマニュアルでは、以下のページに記載があります。また、そのページの末尾からリンクされている EricsBroadcastingDoc には具体例をともなう説明があります。
Broadcasting — NumPy v1.9 Manual
このブロードキャスティングを PHP で実装してみました。コードは以下のとおりです。関数名の bsxfun は Binary Singleton eXpansion FUNction の略です。ブロードキャスティング相当の機能は、MATLAB では bsxfun 関数で提供されており、今回はその関数名を使いました。
<?php function bsxfun($f, $a, $b) { $da = depth($a); $db = depth($b); $a = newAxis($a, $db - $da); $b = newAxis($b, $da - $db); return bsxfunRec($f, $a, $b); } function depth($a) { return is_array($a) ? depth($a[0]) + 1 : 0; } function newAxis($a, $n) { for ($i = 0; $i < $n; ++$i) { $a = array($a); } return $a; } function bsxfunRec($f, $a, $b) { if (!is_array($a)) { return $f($a, $b); } $a = expandIfSingleton($a, count($b)); $b = expandIfSingleton($b, count($a)); return array_map( function ($a, $b) use ($f) { return bsxfunRec($f, $a, $b); }, $a, $b); } function expandIfSingleton($a, $n) { if (count($a) == 1) { return array_fill(0, $n, $a[0]); } return $a; }
簡単なテストプログラムでの実行例は次のとおりです。$a は 1 行 3 列の行列、$b は 3 行 1 列の行列を表しています。これらの配列に対して bsxfun 関数で要素ごとの積を求めると、3 行 3 列の行列が得られます。
$ cat bsxfun_test.php <?php require 'bsxfun.php'; $a = [1, 2, 3]; $b = [[1], [2], [3]]; $c = bsxfun(function ($a, $b) { return $a * $b; }, $a, $b); echo json_encode($c) . "\n"; $ php bsxfun_test.php [[1,2,3],[2,4,6],[3,6,9]]
実装した bsxfun 関数は、まず、引数の二つの配列の次元数をそれぞれ計算し (depth 関数)、次元数が異なる場合には、次元数の小さな側を配列で包んでいくことで、両者の次元数を揃えます (newAxis 関数)。この処理は、NumPy のマニュアルの以下の記述に対応します。この例では、一次元配列の Scale を 1x1x3 の三次元配列に変換してから後続の処理を行うようにします*1。
Arrays do not need to have the same number of dimensions. For example, if you have a 256x256x3 array of RGB values, and you want to scale each color in the image by a different value, you can multiply the image by a one-dimensional array with 3 values. Lining up the sizes of the trailing axes of these arrays according to the broadcast rules, shows that they are compatible:
Image (3d array): 256 x 256 x 3 Scale (1d array): 3 Result (3d array): 256 x 256 x 3
bsxfunRec 関数では、要素数が 1 の次元を拡張しながら (expandIfSingleton 関数)*2、再帰呼び出しによって各要素に $f を適用します。この部分は、マニュアルの以下の記述に対応します。
When either of the dimensions compared is one, the other is used. In other words, dimensions with size 1 are stretched or “copied” to match the other.
In the following example, both the A and B arrays have axes with length one that are expanded to a larger size during the broadcast operation:A (4d array): 8 x 1 x 6 x 1 B (3d array): 7 x 1 x 5 Result (4d array): 8 x 7 x 6 x 5
なお、今回の実装では、入力の配列がブロードキャスティングの条件を満たさない場合を考慮していません。これは二つの場合があります。一つは、両者の要素数がどちらも 2 以上で、かつ、異なっている場合です。NumPy のブロードキャスティングでは、そのような場合にはエラーになると定められています。もう一つは、ジャグ配列が与えられた場合です。NumPy では、矩形配列の場合に限りブロードキャスティングが適用されるようです*3。