Interpret neural network like the linear regression equation such as how much will Y change if we change X1 and keep the other variables fixed In linear regression, assume we get the following equation : Y = 0.8X1+1.9X2+2.4X3+4X4.
We can interpret the linear equation: Keep the other predictors fix, one unit change in X1 will increase 0.8 units in Y in average.
Now, I use the neuralnet package in R and built a neural network model. Then, make a plot like the following picture.

From my perspective, all variables start with W such as W11,W21,W31 are the weight that previous variable or neural. The blue line and letters starts with b such as b11,b21,b31 are the bias weight.
Ignore the other coefficients I did not notice, just focus on the X1, can I write the equation of this neural network plot like:
Y = W31*W21*W11*X1 + W31*W21*W12*X1 + W31*W21*W13*X1 + W32*W21*W11*X1 + W32*W21*W12*X1 + W32*W21*W13*X1 + W32*W22*W11*X1 + W32*W22*W12*X1 + W32*W22*W13*X1
It looks quite ugly, but what I am curious is can I interpret the formula like the previous linear model ? Such as how Y changes if X1 increases one unit.
Or are there any other ways to interpret neural network plot (neuralnet) like this?

This is a nice question, that touches on some interesting points in the history of neural networks (which I can only briefly mention here).
First, what you say is absolutely right if and only if the nodes in your network have linear activation functions, e.g. $b_{1,1} = x_1w_1 + x_2w_2 + x_3w_3$.
However, as explained here, if your network only has linear activation functions, you can simplify the equations so that your multilayer network can be replaced with a simple network with no hidden layers where the input nodes connect directly to the outputs, and this simple network is basically just a linear regression model, $y = x_1b_1 + x_2b_2 + \dots + x_nb_n$.
This is why most useful neural networks use non-linear activation functions, and adding making the model non-linear means the equation you want can no longer be calculated: the effect of changing the value of x1 will depend on the values of all the other inputs.
What you can do, however, is compute the average effect of changing x1 in your training data quite easily: by modifying the values of x1 in your data and calculating the average change in y.
