A Method for Calculating the Derivative of Activation Functions Based on Piecewise Linear Approximation

Abstract

Nonlinear functions are widely used as activation functions in artificial neural networks, which have a great impact on the fitting ability of artificial neural networks. Due to the complexity of the activation function, the computation of the activation function and its derivative requires a lot of computing resources and time during training. In order to improve the computational efficiency of the derivatives of the activation function in the back-propagation of artificial neural networks, this paper proposes a method based on piecewise linear approximation method to calculate the derivative of the activation function. This method is hardware-friendly and universal, it can efficiently compute various nonlinear activation functions in the field of neural network hardware accelerators. In this paper, we use least squares to improve a piecewise linear approximation calculation method that can control the absolute error and get less number of segments or smaller average error, which means fewer hardware resources are required. We use this method to perform a segmented linear approximation to the original or derivative function of the activation function. Both types of activation functions are substituted into a multilayer perceptron for binary classification experiments to verify the effectiveness of the proposed method. Experimental results show that the same or even slightly higher classification accuracy can be achieved by using this method, and the computation time of the back-propagation is reduced by 4–6% compared to the direct calculation of the derivative directly from the function expression using the operator encapsulated in PyTorch. This shows that the proposed method provides an efficient solution of nonlinear activation functions for hardware acceleration of neural networks.

1. Introduction

In recent years, various artificial neural networks (ANNs) accelerators have emerged, including hardware accelerations for convolutional neural networks [1,2,3], and transformer [4,5,6], which all involve the implementation of activation functions, however, these accelerators are basically oriented to the inference only, and there are only a few accelerators for training tasks, such as [7,8]. Above work either uses simple activation functions, such as ReLu, whose derivatives is also simple, or skip the calculation of the derivative of the activation function. The calculation of the derivative of the activation function is often one of the difficult points in the training task of ANNs. Most of the activation functions are complex transcendental functions, such as sigmoid, swish [9], and softplus [10], which provide nonlinear fitting capability for ANNs, and their derivative calculations involve not only derivative operations, but also complex transcendental function calculations. Therefore, it is necessary to optimize the calculation of the derivative of the activation function.

Similar to the way many accelerators implement the activation function, this paper uses a piecewise linear(PWL) approximation method, which requires only one multiply-add to calculate the derivative of the activation function. Activation function can be divided into two types: (1) the derivative function can be represented by the original function, such as sigmoid: $y=1/\left(1+exp\left(-x\right)\right)$, whose derivative function $y^{\prime }=y\left(1-y\right)$; (2) the derivative function can’t be expressed by the original function, such as softplus: $y=log \left(1+exp\left(x\right)\right)$, swish: $y=x/\left(1+exp\left(-x\right)\right)$. For the first type of activation function, the original function value is first calculated using the PWL approximation method, and then the derivative value is calculated indirectly, meanwhile, analyze the transfer error, while for another type of activation function, the derivative function is fitted directly using a PWL approximation method. Then, we add the function obtained by approximation to the neural network in the form of a custom activation function and complete the training test task. In the following, the background knowledge of neural network training and the development of PWL approximation algorithm are briefly introduced.

The training of ANNs is done by adjusting the parameters of the hidden and output layers so that the results computed by the network are as close to the real ones as possible. The training process consists of two parts: forward propagation and backward propagation [11]. In the forward propagation process, the training data is computed through the weights, bias, and activation functions to obtain the hidden layer, and the hidden layer gets the next hidden layer through the weights, bias, and activation functions of the next level, and the output vector, usually the classification result, is finally obtained after layer-by-layer iteration.

The basic principle of back-propagation is that the loss function is first calculated based on the forward propagation output and the label, and then some optimization methods such as gradient descent are used to calculate the bias derivative of the loss function for each weight and bias by the chain rule, the effect of the weight or bias on the loss function, and finally the weights and biases are updated. According to [12], the process of calculating the bias derivatives by the chain rule is actually done by constructing the Jacobi matrix of each layer and then calculating the vector Jacobi product (VJP), where the activation function layer is a point-wise layer and its Jacobi matrix is a diagonal matrix, in practice, it is more efficient to calculate the VJP directly.

PWL is popular in the field of transcendental function approximation calculations due to its computational efficiency and memory friendliness. Many PWL algorithms have been proposed, and the core idea is the trade-off between the number of segments and accuracy, since more segments means more storage space is required. Frenzen et al. [13] studied the number of segments required for approximation calculations of various commonly used functions. The development of segmentation algorithms can be broadly summarized as follows: uniform segmentation, non-uniform segmentation, and adaptive segmentation. Initially, PWL basically used uniform segmentation, and [14] implement the approximate computation of $log$ functions in hardware based on uniform segmentation. Obviously, uniform segmentation has great limitations, and the error gap is large for different functions and different segments. So non-uniform segmentation is proposed, and [15] propose a non-uniform segmentation approach to approximate the fit of the $log$ function, and it divides the whole interval of independent variables into 15 segments and more segments near the zero point to improve the computational accuracy. This approach is still not accurate enough and cannot be extended to other functions.

Recently, the segmentation algorithm gradually tends to a general way with controlled error and fewer segments, [16] propose a method to determine the segmentation points based on the second-order derivatives, it holds that the second-order derivatives reflect the degree of concavity of the function curve, and more segments should be divided where the absolute value of the second-order derivatives is large, this work can make the $L_2$ error below ${10}^{-4}$ with 64 segments, however, the obvious drawback is that it needs to calculate the $2/5$th power of the second order derivative and introduces an integration operation, which makes it difficult to calculate segmentation points. In [17], S.R. Chiluveru et al. use the least squares method to approximate the transcendental function and iterate over the input interval to find the interval that satisfies the error requirement, however, the method is proposed for continuous intervals and needs improvement in hardware implementation. A general PWL method with controllable absolute error is proposed in [18], which determines the slope and intercept of the fitted line by the starting and ending points of the subinterval, and moves the line vertically to control the maximum absolute error. To the best of our knowledge, it is state of the art. On its basis, we use the least squares method, instead of simply using the subinterval start and end points, to determine the slope and intercept of the fitted straight line, which means that all points within the subinterval are taken into account. Experiments show that the improved method, in specific cases, yields a smaller average error or a smaller number of segments. We call the method PWLMMAE (Piecewise Linear Minimize Maximum Absolute Error) and will describe it in detail in the next section.

To summarize, this paper focuses on the computation of the derivative values of the activation function for back-propagation in ANNs training tasks and makes the following contributions.

• We use least squares to improve a general, error-controlled PWL approximation method to obtain fewer segments or smaller average errors, and then extend it to the calculation of various activation functions and their derivatives.
• We evaluate the applicability of the method in neural networks in terms of convergence speed, generalization ability, and hardware overhead.
• We replace the derivative calculation in the neural network training task with proposed method and verified its effectiveness on a three-layer perceptron: our method reduced the backpropagation computation time by 4–6% with the same or even slightly higher classification accuracy.

2. Methods

In this section we first describe the flow of the PWLMMAE algorithm and then apply it to two typical activation functions and analyze their errors.

2.1. PWLMMAE

Similar to [19,20,21], the core idea of the algorithm is to determine the subinterval straight line by least squares, then calculate the maximum absolute error between the line and the real curve, and find the maximum absolute error less than the predetermined error through continuous iteration, the steps of the algorithm are as follows.

2.1.1. Input Range Discretization

Considering the hardware implementation, the input range should be discretized. For a given input interval [M, N], the input $x$ should be defined as a vector

$$x=x\left(1:NUM\right)=M,M+\frac{1}{2^Q},M+\frac{1}{2^Q},\dots ,N$$

(1)

where Q is the number of fractional bits setting in hardware and NUM is the length of the vector.

2.1.2. Minimization of MAE for a Given Width of Subinterval

The slope b and intercept a of the subinterval approximation line are first calculated using the least squares method by Equations (2) and (3), and we can use$h^{\prime }\left(x\right)$ to represent the approximation line.

$$bn+a{\displaystyle \sum }_{i=j}^k{x}_i={\displaystyle \sum }_{i=j}^{k}f\left(x_i\right)$$

(2)

$$b{\displaystyle \sum }_{i=j}^k{x}_i+a{\displaystyle \sum }_{i=j}^k{x}_i^2={\displaystyle \sum }_{i=j}^k{x}_{i}f\left(x_i\right)$$

(3)

where n is the number of discrete points in the subinterval, $x\in \left(j:k\right)$, $1\le j\le k\le NUM$, then, the objective function is denoted by $f\left(x\right)$, so the error vector can be expressed as Equation (4)

$$\delta =f\left(x\left(j:k\right)\right)-h^{\prime }\left(x\left(j:k\right)\right)$$

(4)

The corresponding $MAE$ can also be calculated as

$$MAE=\left\{ \left|max\left(\delta \right)\right|,\left|\min \left(\delta \right)\right| \right\}$$

(5)

2.1.3. Segmentation Points

To obtain the maximum segmentation interval, we determine the segmentation points from right to left. Initially, we set $START=x\left(1\right)$, $END=x\left(NUM\right)$, then perform the PWL method on $x\left(START:END\right)$ to calculate MAE by Equation (4). if $MAE<E_c$, where $E_c$ is a predefined error, approximation succeeds; otherwise $END=END-1$ and repeat above step. Once approximation succeeds, we will update $START$ and $END$ to find the next subinterval, where MAE < E_c. The values of $START$ and $END$ record the segmentation points. This process is shown in Figure 1, and the flow chart of the whole algorithm is as Figure 2.

2.2. PWL Approximation to Derivatives

We choose two typical activation functions, one is sigmoid and the other is softplus, and approximate their derivatives separately setting $E_c=0.001$.

For sigmoid, its function expression and derivative expression are as follows.

$$y=\frac{1}{1+e^{-x}}$$

(6)

$$y^{\prime }=y\left(1-y\right)$$

(7)

We use PWLMMAE to fit its original function and then calculate its derivative by Equation (6). Here, we need to consider the problem of error transmission. When fitting the original function, we can control the absolute error below $E_c$, and this error will be transmitted to the derivative with Equation (7). According to the error propagation law, the error of derivative can be expressed as

$$\delta y^{\prime }=\left(1-2y\right)\delta y$$

(8)

In Equation (8), $\delta y$ is the original error that is smaller than $E_c$, and $y\in \left(0,1\right)$, so that we can find that this calculation does not enlarge the absolute error of the derivative, $\delta y^{\prime }$. It’s also smaller than $E_c$. The same problem will happen to the tanh function, the difference is that the error of derivative of tanh will be enlarged to twice of the original error, however its effects is marginal, and we can still control the error of its derivative within a predefined value. The reason why we do not use PWLMMAE to fit the derivatives directly is that the above approach can simplify both the calculation of the original function in forward propagation and the derivatives in backward propagation with controlled errors and no increase in hardware area.

For softplus, its function expression and derivative expression are as follows.

$$y=log\left(1+e^x\right)$$

(9)

$$y^{\prime }=\frac{e^x}{1+e^x}$$

(10)

We fit its derivative function directly using PWLMMAE, so there is no need to consider the error transfer. Figure 3 shows the approximation of the derivative of softplus and sigmoid. We can find that the max absolute error in the derivative of sigmoid is expanded by a factor of two.

3. Experiment

In this section, we first perform comparison experiments of PWL methods. Then, we perform binary classification experiments on a three-layer perceptron, in which we evaluate the suitability of the PWLMMAE algorithm for implementing derivatives of the activation function in three ways. Finally, we compare the classification accuracy and computational speed of our method with the operator encapsulated in PyTorch.

3.1. Experimental Setup

In the first part, we approximate the two typical activation functions, sigmoid and tanh, and compare with [18]. In the binary classification experiments, we use the ionosphere dataset [22]. We replaced the activation layers in the network with our own implementation of sigmoid and softplus, respectively. For comparison, we implemented two versions of each activation function, one version using the PWLMMAE and the other version using operators such as torch.exp, torch.log in PyTorch. To determine the interval, we counted the input distribution of the activation layers in the network several times, as shown in the Figure 4, which allows us to determine the input range as (−15, 15), and for the very few inputs that are not in that range, we use the nearest segment linear approximation when calculating the derivatives.

3.2. Segment Approximation Comparison

We conduct comparison experiments for the sigmoid function and tanh function at different maximum absolute errors, respectively. The experimental settings and results are shown in the

Table 1. Segment Approximation Comparison.

Function	Input_Range	Q	${\mathit{E}}_{\mathit{c}}$	Method	${\mathit{E}}_{\mathit{a}}$	NO. of Segment
sigmoid	(−15,15)	4	0.001	[18]	4.65 × 10−4	21
				This	3.96 × 10−4	20
			0.0005	[18]	2.61 × 10−4	29
				This	2.16 × 10−4	29
tanh	(−4,4)	4	0.001	[18]	5.00 × 10−4	30
				This	4.94 × 10−4	31
			0.0005	[18]	2.39 × 10−4	52
				This	2.61 × 10−4	46

below. The input range of the sigmoid function is set to (−15, 15), while the input range of the tanh function is set to (−4, 4), beyond which the values of the two functions converge separately. The fraction bit width Q is set to 4 bits. The preset maximum absolute errors are set for two sets, $E_c=0.001$ and $E_c=0.0005$. As can be seen from the Table 1 for the sigmoid function, our method is able to obtain a smaller $E_a$ (average error) with a reduced number of segments when $E_c=0.001$, and for the tanh function, a significant reduction in the number of segments is achieved when $E_c=0.0005$.

3.3. Evaluate PWLMMAE

According to [23], there are three criteria for selecting the implementation of the derivative of activation function: (1) speed of convergence, (2) capability of generalization, and (3) hardware overhead and speed. In the following, we will evaluate the suitability of the PWLMMAE algorithm for implementing the derivative of the activation function in terms of these three aspects.

3.3.1. Speed of Convergence

As can be seen from Figure 5, the loss decreases rapidly in the first 50 epoch. After the same number of epoch, the network using the PWLMMAE method has a smaller loss value. In the last 50 epochs, loss decreases slowly and tends to converge. The sigmoid, softplus functions implemented using the PLWMMAE algorithm converge faster than sigmoid, softplus function encapsulated in PyTorch.

3.3.2. Capability of Generalization

We use the confusion matrix and the F-score to evaluate the model. The

Table 2. Evaluate sigmoid.

Activation Function	Macro-P	Macro-R	Macro-F1
sigmoid-PyTorch	0.90195	0.82159	0.85989
sigmoid-PWLMMAE	0.90438	0.81338	0.85647

and

Table 3. Evaluate softplus.

Activation Function	Macro-P	Macro-R	Macro-F1
softplus-PyTorch	0.92694	0.79302	0.85477
softplus-PWLMMAE	0.91681	0.80705	0.85844

below shows the macro precision macro-p, macro recall macro-R and macro-F1 calculated from the confusion matrix of 100 training, testing. The model using this method has a very small drop in macro-F1 compared to the model using the PyTorch activation layer, which means that using our method does not significantly reduce the generalization ability of the model.

3.3.3. Hardware Overhead and Speed

It is obvious that only one multi-add operation is required to complete a derivative calculation using PWLMMAE. In [18], by designing the indexing circuit properly and using pipelining, it just takes one clock cycle to complete one derivative calculation and the hardware area is significantly better than other PWL methods. In addition, we study the accuracy of the input interval discretization, i.e., the bit width of the fractional part of the hardware implementation. In the hardware implementation of neural networks, aggressive reduction of the data bit width can reduce the storage and operation overhead to a great extent [24], and for PWL approximation, it can also reduce the number of segments and reduce the storage overhead. As shown in the

Table 4. Effect of bit width of data on classification accuracy.

Function	Q(Bits)	${\mathit{E}}_{\mathit{c}}$	Acc Average	NO. of Segment
Sigmoid	2	0.25	87.27	2
	3	0.125	88.93	3
	4	0.0625	89.82	3
Softplus	2	0.25	90.10	2
	3	0.125	90.01	3
	4	0.0625	90.19	3

below, when the bit width is reduced, the corresponding predefined maximum error $E_c$ changes to accommodate the max precision that can be represented by the bit width. For the sigmoid function, the classification accuracy of the neural network decreases when the bit width of the fractional part is reduced, while for the softplus function, this does not occur. We analyze the reason is that in the custom implementation of both activation functions, the sigmoid original function is approximated, while the derivative function of the softplus function is approximated, which means that the value of its original function is still exact.

3.4. Performance Comparisons

We train the network 100 times and count the average backward time during training, while verifying the results of each training on the test dataset and recording the Distribution of inputs for each activation layer classification accuracy. The Backward time refers to the time used to take the reverse derivative from the loss function and update the weights of each layer, which is called back-propagation. From the

Table 5. Evaluate sigmoid.

Activation Function	Accuracy	Backward Time(s)
Sigmoid-PyTorch	90.00	2.031
Sigmoid-PWLMMAE	89.82	1.935
Softplus-PyTorch	89.91	2.139
Softplus-PWLMMAE	90.19	2.005

, we can see that our method can reduce the backward time by more than 4% for the sigmoid function and about 6% for the softplus function, while guaranteeing almost no impact on the classification accuracy.

4. Conclusions

This paper focuses on the calculation of nonlinear activation functions in ANNs. We propose a generalized, error-controlled PWL approximation method PWLMMAE using least squares, which is capable of obtaining smaller average approximation errors with fewer segments. This method can calculate the approximate calculation of any nonlinear function quickly with less hardware resources. We explore the possibility of applying this method to neural network hardware acceleration and use it to calculate the derivative values of the activation functions in the neural network training task. We first evaluate the applicability of the method in ANNs in terms of convergence speed, generalization capability, hardware area and speed. Experimental results show that the activation function implemented using the method is almost indistinguishable from the activation function encapsulated in PyTorch and is capable of performing well. Finally, we compare the accuracy, back-propagation computation time of the activation function implemented by this method and PyTorch operator on the test dataset, and our method reduces the back-propagation computation time by 4–6% with the same or even slightly higher classification accuracy. Based on the above experiments, it can be seen that the proposed approximate calculation method is suitable for the hardware acceleration of ANNs, including reasoning and training. However, running the algorithm at the software level is limited by the indexing time when calculating the derivatives, and we have only experimented on a three-layer perceptron, future work will see us apply the method to larger-scale networks and complete the hardware implementation.