I am a cognitive scientist, and am modelling a particular kind of neural network called a "masking field", where the change in activation of a particular node in the masking field is:
$${dy_i(t)\over dt} = -Ay_i(t) + [B - y_i(t)]P_i(t) - [C + y_i(t)]Q_i(t)$$
Where $A$, $B$ and $C$ are constants, $y_i(t)$ is the activation of node $i$, $P_i(t)$ is the net excitation to this node from other nodes, and $Q_i(t)$ is the net inhibition from other nodes. Excitation/inhibition from other nodes is dependent on their activation levels, which also change with time.
It is straightforward enough to model this iteratively in computer using Euler's method:
$$y_{in+1} = y_{in} + h\{-Ay_{in} + [B - y_{in}]P_{in} - [C + y_{in}]Q_{in}\}$$
where $h$ is the step size, and all $y_{i0}$ are known. For each step, I calculate all the $P_{in}$ and $Q_{in}$ for all nodes using node activations at timestep $n$, then calculate the activation values of all nodes for step $n+1$. The next step I can use all the updated activations to calculate all of the $P_{in+1}$ and $Q_{in+1}$, before calculating the $n+2$ activations, etc.
However, I wanted to see if using the Runge-Kutta fourth order method would give better results. My question is:
Q. When calculating the K1, K2, K3 and K4 values for a particular node $i$ using the Runge-Kutta method, do I recalculate $P_{in}$ and $Q_{in}$ when calculating each of the four K values, since the activations of all the other nodes are also changing? That is, I would calculate K1 for every node, and use all these K1 values to then calculate K2 values for every node, then use all these to calculate the K3s, and the K4s, before finally calculating the n+1 activations? Or should I treat $P_{in}$ and $Q_{in}$ as constant when calculating K1-4 for node $i$, so that I am only considering the change in activation of one node at a time?
Thankyou!
If you were to treat $P$ and $Q$ as constant during the step, then you would be better served by using the exact solution of $$ \dot y+ay=b,\qquad a=A+P+Q, \;b=BP-CQ $$ where the solution is $$ y(t+h)=b+e^{-ah}(y(t)-b) $$ which can be then used instead of the Euler step.
For the full Runge-Kutta method you got it right, you have to re-compute the $P$ and $Q$ in every stage. However, you can use much larger step-sizes than with the Euler method for the same accuracy. Improving the accuracy by a factor of $16=2^4$ requires 16 times the steps with the Euler method, but only a doubling of steps and thus function evaluations with RK4.