1function Y = dnn (W, bias, Y0)
2%GRB.DNN Sparse deep neural network in GraphBLAS.
3% Performs ReLU inference using input feature vector(s) Y0, DNN weights W,
4% and bias vectors.  The input features are in a matrix Y0 of size
5% nfeatures-by- nneurons.  The DNN weights W is a cell array with W{k}
6% being the kth layer of the DNN, so that the number of layers is nlayers =
7% length (W).  W{k} is a matrix of size nneurons-by-nneurons.  The bias
8% variable is a cell array of length nlayers.  Each bias{k} is a diagonal
9% matrix of size nneurons-by-nneurons, which gives the bias values of each
10% neuron in the kth layer.
11%
12% Usage:
13%
14%   Y = GrB.dnn (W, bias, Y0) ;
15%
16% The matrices can be stored by row or by column, but GrB.format ('by row')
17% is somewhat faster.  For the 2019 GraphChallenge, all matrices can be
18% 'single', and the same results are obtained.
19%
20% In the MATLAB reference implementation, the bias{k} is a row vector of
21% size 1-by-nneurons.  The MATLAB reference inputs can be converted to
22% GraphBLAS matrices with the following code:
23%
24%   d = struct ('format', 'by row') ;
25%   n = size (Y0, 2) ;
26%   Y0 = GrB (Y0, 'single', 'by row') ;
27%   for k=1:length(W)
28%       W {k} = GrB (W {k}, 'single', 'by row') ;
29%       bias {k} = GrB.build (1:n, 1:n, bias {k}, n, n, '+', 'single', d) ;
30%   end
31%
32% All of the above conversion is optional, except for bias {k} since it is
33% changed from a row vector to a diagonal matrix.
34%
35% See also dnn_matlab, dnn_mat2gb.
36
37% SuiteSparse:GraphBLAS, Timothy A. Davis, (c) 2017-2021, All Rights Reserved.
38% SPDX-License-Identifier: GPL-3.0-or-later
39
40% NOTE: this is a high-level algorithm that uses GrB objects.
41
42[f,~] = GrB.format (Y0) ;
43desc.format = '' ;
44if (isequal (f, 'by row'))
45    % hypersparse-by-row is fastest, since entire rows drop out of Y
46    desc.format = 'hyper by row' ;
47end
48tol = single (32) ;
49
50Y = Y0 ;
51for k = 1:length(W)
52    % Propagate through layer, apply bias, and threshold negative values.
53    Y = GrB.mxm (Y, '+.*', W {k}, desc) ;
54    Y = GrB.select (GrB.mxm (Y, '+.+', bias {k}, desc), '>0', desc) ;
55    M = Y > tol ;
56    if (nnz (M) > 0)
57        % Y (M) = tol ;
58        Y = GrB.subassign (Y, M, tol, desc) ;
59    end
60end
61
62