Answers.Body: <p>First, a forward pass is done to obtain predictions of <span class="math-container">$z$</span>, at every <span class="math-container">$t$</span>. Then the adjoint state is run backward in time for every <span class="math-container">$t$</span>. Which gives the learning impulse. So an initial run is done to obtain values of the dynamical system at all time points and the last of these is the initial point for the backward pass.</p>
  <p>This picture depicts it.
  <a href="https://i.sstatic.net/qvtfy.png" rel="nofollow noreferrer"><img src="https://i.sstatic.net/qvtfy.png" alt="enter image description here" /></a></p>
  <p>To be more specific with regards to the question in your initial paragraph the input of <span class="math-container">$t_0$</span> and <span class="math-container">$t_1$</span> into <span class="math-container">$ODE\_Solve$</span> are time points at which the network is evaluated. In this specific case, the first point is <span class="math-container">$t_0$</span> which should always be the case and <span class="math-container">$t_1$</span> is the last point. But you can also insert more time points for example <span class="math-container">$N$</span> time points than you input the sequence <span class="math-container">$t_0, t_1, t_2, ..., t_N$</span> into <span class="math-container">$ODE\_Solve$</span>. In that case, <span class="math-container">$t_N$</span> is the last time point.</p>
  <p>This picture also shows a lot. But in this picture, a RNN is used to generate <span class="math-container">$z_{t_0}$</span>. The RNN encoder is given the time points in reversed order: <span class="math-container">$(t_N, t_{N-1}, t_{N-2}, ..., t_0)$</span>.
  <a href="https://i.sstatic.net/YnJnC.png" rel="nofollow noreferrer"><img src="https://i.sstatic.net/YnJnC.png" alt="enter image description here" /></a></p>
  <p>This is probably not enough to satisfy your question, hopefully, I can get more specific with respect to the math later.</p>
  
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2020-10-16T20:57:05.587Z
Answers.LastActivityDate: 2022-07-07T08:13:48.823Z
Answers.LastEditDate: 2022-07-07T08:13:48.823Z
Answers.UserDisplayName: Quibus
Body: <p>I'm reading the paper <a href="https://arxiv.org/pdf/1806.07366.pdf" rel="nofollow noreferrer">Neural Ordinary Differential Equations</a> and I have a simple question about adjoint method. When we train NODE, it uses a blackbox ODESolver to compute gradients through model parameters, hidden states, and time. It uses another quantity <span class="math-container">$\mathbf{a}(t) = \partial L / \partial \mathbf{z}(t)$</span> called <em>adjoint</em>, which also satisfies another ODE. As I understand, the authors build a single ODE that computes all the gradients <span class="math-container">$\partial L / \partial \mathbf{z}(t_{0})$</span> and <span class="math-container">$\partial L / \partial \theta$</span> by solving that single ODE. However, I can't understand how do we know the value <span class="math-container">$\partial L / \partial \mathbf{z}(t_1)$</span> which corresponds to the initial condition for the ODE corresponds to the adjoint. I'm using <a href="https://msurtsukov.github.io/Neural-ODE/" rel="nofollow noreferrer">this</a> tutorial as a reference, and it defines custom forward and backward methods for solving ODE. However, for the backward computation (especially <code>ODEAdjoint</code> class in the tutorial) we need to pass <span class="math-container">$\partial L / \partial \mathbf{z}$</span> for backpropagation, and this enables us to compute <span class="math-container">$\partial L / \partial \mathbf{z}(t_i)$</span> from <span class="math-container">$\partial L / \partial \mathbf{z}(t_{i+1})$</span>, but we still need to know the adjoint value <span class="math-container">$\partial L / \partial \mathbf{z}(t_N)$</span>. I do not understand well about how pytorch's <code>autograd</code> package works, and this seems to be a barrier to understand this. Could anyone explain how it operates, and where <span class="math-container">$\partial L / \partial \mathbf{z}(t_1)$</span> (or <span class="math-container">$\partial L / \partial \mathbf{z}(t_N)$</span> if this is more comfortable) comes from? Thanks in advance.</p>
  <hr />
  <p>Here's my <em>guess</em> for the initial adjoint from simple example. Let <span class="math-container">$d\mathbf{z}/dt = Az$</span> be a 2-dim linear ODE with given <span class="math-container">$A \in \mathbb{R}^{2\times 2}$</span>. If we use Euler's method as a ODE solver, then the estimate for <span class="math-container">$z(t_1)$</span> is explicitly given as <span class="math-container">$$\hat{\mathbf{z}}(t_1) = \mathrm{ODESolve}(\mathbf{z}(t_0), f, t_0, t_1, \theta))= \left(I + \frac{t_1 - t_0}{N}A\right)^{N} \mathbf{z}(t_0) $$</span> where <span class="math-container">$N$</span> is the number of steps for Euler's method (so that <span class="math-container">$h = (t_1 - t_0) /N$</span> is the step size). If we use MSE loss for training, then the loss will be
  <span class="math-container">$$
  L(\mathbf{z}(t_1)) = \Bigl|\Bigl| \mathbf{z}_1 - \left(I + \frac{t_1 - t_0}{N}A\right)^N\mathbf{z}(t_0)\Bigr|\Bigr|_2^2
  $$</span>
  where <span class="math-container">$\mathbf{z}_1$</span> is the true value at time <span class="math-container">$t_1$</span>, which is <span class="math-container">$\mathbf{z}_1 = e^{A(t_1 - t_0)}\mathbf{z}(t_0)$</span>. Since adjoint <span class="math-container">$\mathbf{a}(t) = \partial L / \partial \mathbf{z}(t)$</span> satisfies <span class="math-container">$$\frac{d\mathbf{a}(t)}{dt} = -\mathbf{a}(t)^{T} \frac{\partial f(\mathbf{z}(t), t, \theta)}{\partial \mathbf{z}} = \mathbf{0},$$</span>
  <span class="math-container">$\mathbf{a}(t)$</span> is constant and we get <span class="math-container">$\mathbf{a}(t_0) = \mathbf{a}(t_1)$</span>. So we do not need to use augmented ODE for computing <span class="math-container">$\mathbf{a}(t)$</span>. However, I still don't know what <span class="math-container">$\mathbf{a}(t_1) = \partial L / \partial \mathbf{z}(t_1)$</span> should be.  If my understanding is correct, since <span class="math-container">$L = ||\mathbf{z}_1 - \mathbf{z}(t_1)||^{2}_{2}$</span>, it seems that the answer might be
  <span class="math-container">$$
  \frac{\partial L}{\partial \mathbf{z}(t_1)} = 2(\mathbf{z}(t_1) - \mathbf{z}_1).
  $$</span>
  However, this doesn't seem to be true: if it is, and if we have multiple datapoints at <span class="math-container">$t_1, t_2, \dots, t_N$</span>, then the loss is
  <span class="math-container">$$
  L = \frac{1}{N} \sum_{i=1}^{N}||\mathbf{z}_i  -\mathbf{z}(t_i)||_{2}^{2}
  $$</span>
  and we may have
  <span class="math-container">$$
  \frac{\partial L}{\partial \mathbf{z}(t_i)} = \frac{2}{N} (\mathbf{z}(t_i) - \mathbf{z}_i),
  $$</span>
  which means that we don't need to solve ODE associated to <span class="math-container">$\mathbf{a}(t)$</span>.</p>
  
Tags: <neural-networks><backpropagation>
Title: Computation of initial adjoint for NODE
UserDisplayName: Seewoo Lee
ThreadText: # Computation of initial adjoint for NODE
  
  (Asked by: username_0 on 2020-10-12)
  
  I'm reading the paper [Neural Ordinary Differential Equations](https://arxiv.org/pdf/1806.07366.pdf) and I have a simple question about adjoint method. When we train NODE, it uses a blackbox ODESolver to compute gradients through model parameters, hidden states, and time. It uses another quantity $\mathbf{a}(t) = \partial L / \partial \mathbf{z}(t)$ called *adjoint*, which also satisfies another ODE. As I understand, the authors build a single ODE that computes all the gradients $\partial L / \partial \mathbf{z}(t\_{0})$ and $\partial L / \partial \theta$ by solving that single ODE. However, I can't understand how do we know the value $\partial L / \partial \mathbf{z}(t\_1)$ which corresponds to the initial condition for the ODE corresponds to the adjoint. I'm using [this](https://msurtsukov.github.io/Neural-ODE/) tutorial as a reference, and it defines custom forward and backward methods for solving ODE. However, for the backward computation (especially `ODEAdjoint` class in the tutorial) we need to pass $\partial L / \partial \mathbf{z}$ for backpropagation, and this enables us to compute $\partial L / \partial \mathbf{z}(t\_i)$ from $\partial L / \partial \mathbf{z}(t\_{i+1})$, but we still need to know the adjoint value $\partial L / \partial \mathbf{z}(t\_N)$. I do not understand well about how pytorch's `autograd` package works, and this seems to be a barrier to understand this. Could anyone explain how it operates, and where $\partial L / \partial \mathbf{z}(t\_1)$ (or $\partial L / \partial \mathbf{z}(t\_N)$ if this is more comfortable) comes from? Thanks in advance.
  
  ---
  
  Here's my *guess* for the initial adjoint from simple example. Let $d\mathbf{z}/dt = Az$ be a 2-dim linear ODE with given $A \in \mathbb{R}^{2\times 2}$. If we use Euler's method as a ODE solver, then the estimate for $z(t\_1)$ is explicitly given as $$\hat{\mathbf{z}}(t\_1) = \mathrm{ODESolve}(\mathbf{z}(t\_0), f, t\_0, t\_1, \theta))= \left(I + \frac{t\_1 - t\_0}{N}A\right)^{N} \mathbf{z}(t\_0) $$ where $N$ is the number of steps for Euler's method (so that $h = (t\_1 - t\_0) /N$ is the step size). If we use MSE loss for training, then the loss will be
  $$
  L(\mathbf{z}(t\_1)) = \Bigl|\Bigl| \mathbf{z}\_1 - \left(I + \frac{t\_1 - t\_0}{N}A\right)^N\mathbf{z}(t\_0)\Bigr|\Bigr|\_2^2
  $$
  where $\mathbf{z}\_1$ is the true value at time $t\_1$, which is $\mathbf{z}\_1 = e^{A(t\_1 - t\_0)}\mathbf{z}(t\_0)$. Since adjoint $\mathbf{a}(t) = \partial L / \partial \mathbf{z}(t)$ satisfies $$\frac{d\mathbf{a}(t)}{dt} = -\mathbf{a}(t)^{T} \frac{\partial f(\mathbf{z}(t), t, \theta)}{\partial \mathbf{z}} = \mathbf{0},$$
  $\mathbf{a}(t)$ is constant and we get $\mathbf{a}(t\_0) = \mathbf{a}(t\_1)$. So we do not need to use augmented ODE for computing $\mathbf{a}(t)$. However, I still don't know what $\mathbf{a}(t\_1) = \partial L / \partial \mathbf{z}(t\_1)$ should be. If my understanding is correct, since $L = ||\mathbf{z}\_1 - \mathbf{z}(t\_1)||^{2}\_{2}$, it seems that the answer might be
  $$
  \frac{\partial L}{\partial \mathbf{z}(t\_1)} = 2(\mathbf{z}(t\_1) - \mathbf{z}\_1).
  $$
  However, this doesn't seem to be true: if it is, and if we have multiple datapoints at $t\_1, t\_2, \dots, t\_N$, then the loss is
  $$
  L = \frac{1}{N} \sum\_{i=1}^{N}||\mathbf{z}\_i -\mathbf{z}(t\_i)||\_{2}^{2}
  $$
  and we may have
  $$
  \frac{\partial L}{\partial \mathbf{z}(t\_i)} = \frac{2}{N} (\mathbf{z}(t\_i) - \mathbf{z}\_i),
  $$
  which means that we don't need to solve ODE associated to $\mathbf{a}(t)$.
  
  **Tags:** neural-networks, backpropagation
  
  
  ---
  
  ## 1 answers
  
  ### Answer by username_1 on 2020-10-16. Score: 1
  
  First, a forward pass is done to obtain predictions of $z$, at every $t$. Then the adjoint state is run backward in time for every $t$. Which gives the learning impulse. So an initial run is done to obtain values of the dynamical system at all time points and the last of these is the initial point for the backward pass.
  
  This picture depicts it.
  [![enter image description here](https://i.sstatic.net/qvtfy.png)](https://i.sstatic.net/qvtfy.png)
  
  To be more specific with regards to the question in your initial paragraph the input of $t\_0$ and $t\_1$ into $ODE\\_Solve$ are time points at which the network is evaluated. In this specific case, the first point is $t\_0$ which should always be the case and $t\_1$ is the last point. But you can also insert more time points for example $N$ time points than you input the sequence $t\_0, t\_1, t\_2, ..., t\_N$ into $ODE\\_Solve$. In that case, $t\_N$ is the last time point.
  
  This picture also shows a lot. But in this picture, a RNN is used to generate $z\_{t\_0}$. The RNN encoder is given the time points in reversed order: $(t\_N, t\_{N-1}, t\_{N-2}, ..., t\_0)$.
  [![enter image description here](https://i.sstatic.net/YnJnC.png)](https://i.sstatic.net/YnJnC.png)
  
  This is probably not enough to satisfy your question, hopefully, I can get more specific with respect to the math later.
  
  
  ---
  
  