# Which optimization method can be used to do the following?

I've the following system of equations for studying information flow in the below graph,

$$\frac{d \phi}{dt} = -M^TDM\phi + \text{noise effects} \hspace{1cm} (1)$$ Here, M is the incidence matrix of the graph

$$\phi$$ is a vector with variables [ A B C D E F].

I've solved the above ode's to obtain the time series data of variables A,B,C,D,E,F.

Using the time-series data obtained from the above step, I'd like to do determine $$\tilde{D}$$ for the following system

$$\frac{d \phi}{dt} = -M^T\tilde{D}M\phi \hspace{1cm} (2)$$

Note: The entries in the diagonal elements of $$\tilde{D}$$ are the edge weights.

I'd like to know if there is a way to solve this as an optimization problem using time-series data as input.

Any suggestions on how to proceed will be really helpful

EDIT: Equation (1) (with noise effects) is solved using prior values of the diagonal matrix, D and the time series profiles of variables in each node are obtained.

What do I want to achieve?

I want to determine a modified D i.e $$\tilde{D}$$ that can generate the same time series profile that was generated while solving equation (1).

EDIT: Based on the below discussion

The cost function that I want to minimize is

$$\int_{0}^{t_{f}} || \phi(t) - \exp{(-M^{T} \tilde{D} M t)} \phi(0) ||^{2} dt$$

to reduce the difference between the solutions of (1) and (2).

Any suggestions on the optimization function and toolbox that can be used to solve the above will be really helpful

EDIT3: I looked at some literature and I think this problem has to be solved using fmincon and it might not be possible to use convex optimization.

I'm still looking for suggestions on how to formulate the constraints for solving the cost function illustrated below.

• It’s not quite clear to me what it is you are trying to achieve. What is the cost function you would like to optimize? It appears you are trying to find a system of ODEs that generates some calculated time-series data. If that’s the case, there are several methods. If you take the Laplace transform of your ODE, the problem is transformed into a rational fitting problem, for which several methods exist. A popular one is called vector-fitting. There are also direct time-domain methods, such as time-domain vector fitting. Mar 15, 2020 at 22:09
• @AmitHochman I'm trying to find a modified D i.e $\tilde{D}$ that can generate the calculated time-series data from equation(1). I 've also made an edit in my post to explain it further. Please let me know if it is still unclear. I'll have a look at the time-domain vector fitting method that you mentioned. Will this method be useful for finding the diagonal entris of the matrix $\tilde{D}$ Mar 16, 2020 at 2:42
• @AloneProgrammer Thanks a lot for the response. Could you please mathematically illustrate a bit on why the second equation will not yield a unique solution for $-M^{T} (\tilde{D}-D) M$ ? Mar 16, 2020 at 4:26
• @AloneProgrammer Thanks. I don't completely understand why 36 unknowns. For example, if we consider a particular instance of time let's say $t = t_{i}$, the product of $-M^{T} (\tilde{D}-D) M \phi(t)$ will give a 6 x 1 vector with 6 unknowns (the diagonal entries of $\tilde{D}$), right? Please let me know if I am wrong Mar 16, 2020 at 4:47
• @Natasha I posted a complete answer for you. Mar 16, 2020 at 4:52

Despite my comment, I think you can find $$\tilde{D}$$ that contains the noise term as well. You have this equation:

$$-M^{T} \tilde{D} M \phi(t) = -M^{T} D M \phi(t) + W(t)$$

Where $$W(t)$$ is the noise term vector. So:

$$-M^{T} (\tilde{D}-D) M \phi(t) = W(t)$$

Take $$\mathcal{D} = \tilde{D} - D$$.

Let's expand this equation:

$$(M \phi(t))_{i} = \sum_{j=1}^{6} M_{ij} \phi_{j}$$

$$(\mathcal{D} M \phi(t))_{i} = \sum_{j=1}^{6} \sum_{k=1}^{6} \mathcal{D}_{ij} M_{ik} \phi_{k}$$

But $$\mathcal{D}$$ is diagonal, which means:

$$\mathcal{D}_{ij} = \mathfrak{D}_{ii}\delta_{ij}$$

Where $$\delta_{ij}$$ is the Kronecker delta. So:

$$(\mathcal{D} M \phi(t))_{i} = \sum_{j=1}^{6} \sum_{k=1}^{6} \mathcal{D}_{ij} M_{ik} \phi_{k} = \sum_{j=1}^{6} \sum_{k=1}^{6} \mathfrak{D}_{ii} \delta_{ij} M_{ik} \phi_{k} = \mathfrak{D}_{ii} \sum_{k=1}^{6} M_{ik}\phi_{k}$$

Finally:

$$(M^{T} \mathcal{D} M \phi(t))_{i} = \sum_{j=1}^{6} M_{ji} \mathfrak{D}_{jj} \sum_{k=1}^{6} M_{jk}\phi_{k} = \sum_{j=1}^{6} \sum_{k=1}^{6} M_{ji}M_{jk} \phi_{k} \mathfrak{D}_{jj}$$

or:

$$-(M^{T} \mathcal{D} M \phi(t))_{i} = (W(t))_{i}$$

so:

$$\sum_{j=1}^{6} \sum_{k=1}^{6} M_{ji}M_{jk} \phi_{k} \mathfrak{D}_{jj} = -W_{i}$$

or in matrix form by defining the vector $$\mathsf{x} = (\mathfrak{D}_{11},\mathfrak{D}_{22},\mathfrak{D}_{33},\mathfrak{D}_{44},\mathfrak{D}_{55},\mathfrak{D}_{66})^{T}$$ and matrix $$\mathsf{A}_{ij} = M_{ji}\sum_{k=1}^{6} M_{jk} \phi_{k}$$:

$$\mathsf{A} \mathsf{x} = -W$$

Note that $$\mathsf{A}$$ and $$W$$ depends on time, so for each time you need to solve this equation to extract $$\mathsf{x}$$. In fact due to the fact that $$\tilde{D}$$ suppose to contain the noise term as well, it must depend on time.

So, you can easily find $$\mathsf{x}$$ by solving this linear equation and it would be a unique solution. The only comment here is that: this approach is somewhat looks like a reverse engineering and for a calculated $$\phi$$ and having that particular $$W(t)$$ noise term, you embed it into the diagonal matrix $$D$$.

Update: Another approach to not using an exact method is compromising the accuracy of $$\phi(t)$$ to find a constant $$\tilde{D}$$. The exact equation is:

$$\phi^{'}(t) = -M^{T} D M \phi(t) + W(t)$$

and the approximate equation is:

$$\hat{\phi}^{'}(t) = -M^{T} \tilde{D} M \hat{\phi}(t)$$

Now we want to minimize this:

$$\mathsf{K} = \int_{0}^{t_{f}} ||\phi(t) - \hat{\phi}(t)||^{2} dt$$

Analytically for $$\hat{\phi}(t)$$:

$$\hat{\phi}(t) = \exp{(-M^{T} \tilde{D} M t)} \hat{\phi}(0)$$

So finally your optimization problem is:

$$\tilde{D}, \hat{\phi}(0) = \text{argmin} \ \mathsf{K}(\tilde{D},\hat{\phi}(0)) = \text{argmin} \ \int_{0}^{t_{f}} ||\phi(t) - \exp{(-M^{T} \tilde{D} M t)} \hat{\phi}(0)||^{2} dt$$

This problem could be solved by CVXPY tool I believe.

• Thanks a lot for the detailed solution.IIUC $\mathsf{A} \mathsf{x} = -W$ is for a given instant in time. Mar 16, 2020 at 5:07
• If you insist on finding a constant $\tilde{D}$, you need to do some compromise with respect to the accuracy of $\phi(t)$. For example, formulate a least square or optimization problem as finding $\mathcal{D}$ as a diagonal matrix that minimizes this function for example: $$\int_{0}^{t_{f}} || M^{T} \mathcal{D} M \phi(t) + W(t) ||^{2} dt$$ where $t_{f}$ is your final time. Still, I don't see what you get here by doing this. Is it possible to elaborate why it is necessary to find $\tilde{D}$ and it must be a constant matrix? Mar 16, 2020 at 5:35
• @Natasha You're welcome! Please consider accepting this answer if it solved your problem. Mar 16, 2020 at 6:23
• Yes, definitely :). If you don't mind, could you please suggest a toolbox that can be used to solve the above optimization problem? I've little experience with this, and I am not sure if fmincon of MATLAB will be a suitable choice. Mar 16, 2020 at 6:28
• Sorry, I just saw that you have already mentioned one. Thanks! Mar 16, 2020 at 6:30