# Difficulties to implement the layer-wise relevance propagation in MATLAB

I'm having serious issues with the implementation of the LRP algorithm for neural networks in MATLAB. The challenge is to implement the equations correctly. I'm trying to implement the deep-Taylor $$\alpha_1\beta_0$$ version of the LRP. I'm testing it on a feed-forward full-connected neural network with one hidden layer trained on the MNIST dataset. I start initializing the relevances of the last layer with the predictions of the network:

R{l} = nn.layers{nn.lastLayer}.output;


Then, i need to perform the following equation through the hidden layers (in my case is just one): $$R_i^l = \sum_j \frac{a_i^l \cdot w_{ij}^{+}}{\sum_k a_j^{l+1} \cdot w_{kj}^{+}} \cdot R_j^{l+1}$$

Thanks to a paper, i've followed some suggestions and correctly implemented that equation as follows:

    for l = nn.numOfLayers-1 : -1 : 1
z = max(nn.layers{l+1}.W, 0) * nn.layers{l}.output;
s = R{l+1} ./ z;
c = (max(nn.layers{l+1}.W, 0))' * s;
R{l} = nn.layers{l}.output .* c;
end


I hope the syntax is clear, "nn" is the neural network struct and other fields should be intuitive to understand. The elements are columns, it means that outputs $$a$$ are $$numberOfNeurons \times 1$$. The weights are matrices such that on the rows there are the neurons of the next layer, and on the columns the neurons of the current layer.

Now, my problem is, how do i compute the relevances for the first layer? (input layer) Because, in the input layer i haven't the $$a$$ values (that are the outputs transformed by the activation function). So i've looked also in a paper and found that the following equation can be used to compute the relevances of the first layer (the layer of interest):

$$R_i^1 = \sum_j \frac{w_{ij}^{2}}{\sum_k w_{kj}^{2}} \cdot R_j^{2}$$

I've tried different ways to implement that in matlab, but the results seems wrong. For instance, i've tried with the following code:

for i = 1 : inputSize
scalar = 0;
for j = 1 : nn.layers{1}.layerSize
norm = 0;
for k = 1 : inputSize
norm = norm + nn.layers{1}.W(j,k)^2;
end
scalar = (nn.layers{1}.W(j,i)^2 * R{1}(j))/norm;
Rfin(i) = Rfin(i) + scalar;
end
end


and then return $$R$$. These are the results i should see:

this is an example of what i get:

Do you know what am i doing wrong? Do you have ideas on how to implement correctly in MATLAB the second equation?

Thank you very much.

• Hi and welcome to this community! Can you please ask this question on Data Science SE or Stack Overflow, and then delete this question from this website, because it is off-topic here. – nbro May 19 at 17:19