Data assimilation in price formation

TL;DR

This paper develops a method to estimate buyer and vendor densities in a nonlinear price formation model using measurements, formulating it as an optimal control problem with stability analysis and numerical validation.

Contribution

It introduces a novel approach to estimate densities in price formation models through optimal control, including stability estimates and numerical experiments.

Findings

Stability estimates for controls and densities are established.

Numerical experiments support the theoretical results.

The method effectively estimates densities from measurement data.

Abstract

We consider the problem of estimating the density of buyers and vendors in a nonlinear parabolic price formation model using measurements of the price and the transaction rate. Our approach is based on a work by Puel et al., see \cite{Puel2002}, and results in a optimal control problem. We analyse this problems and provide stability estimates for the controls as well as the unknown density in the presence of measurement errors. Our analytic findings are supported with numerical experiments.

Figures

Equations

$${\partial }_{t}f-\frac{{\sigma }^2}{2}{\partial }_{xx}f=\Lambda (t)({\delta }_{p(t)-a}-{\delta }_{p(t)+a}),x\in \Omega ,t>0,$$

$$\Lambda (t)=-\frac{{\sigma }^2}{2}{\partial }_{x}f(p(t),t),f(p(t),t)=0,$$

$$f(x,0)=f_0(x),p(0)=p_0.$$

$$f_0(p_0)=0,f_0(x)>0\text{ for }x<p_0\text{ and }{f}_0(x)<0\text{ for }x>p_0,\text{ a.e. in }\Omega$$

$${\partial }_{x}f=0\text{ on }\partial \Omega$$

$$F(x,t)$$

$${\partial }_{t}F(x,t)-{\partial }_{xx}F(x,t)$$

$$F(x,0)$$

$$F_0(x)$$

$$F_0(x,t)$$

$$k_L=(p_0+L)/a\text{ and }{k}_R=(L-p_0)/a.$$

$$f(x,t)=F(x,t)-F^{+}(x+a,t)+F^{-}(x-a,t).$$

$${\partial }_{x}F(-L,t)$$

$${\partial }_{x}F(L,t)$$

$$p_{\infty }=\frac{2M^{l}L-a(M^l-M^r)}{2(M_l+M_r)}-\frac{L}{2}$$

$$p^{\prime }(t)=\frac{{\partial }_{t}F(p(t),t)}{\Lambda (t)}=\frac{{\partial }_{xx}F(p(t),t)}{\Lambda (t)},$$

$$\underset{t\in [\epsilon ,T]}{\sup }{p}^{\prime }(t)\le \frac{\Vert {\partial }_{xx}F{\Vert }_{C([-L,L]\times [\epsilon ,T])}}{\underline{\Lambda }},$$

$${\Omega }_{◃}=[-L,p(t)],{\Omega }_{▹}=[p(t),L],\text{ and }\Omega =[0,1],$$

$$Q_{◃}={\Omega }_{◃}\times [\epsilon ,T],Q_{▹}={\Omega }_{▹}\times [\epsilon ,T],\text{ and }Q=\Omega \times [\epsilon ,T].$$

$${\varphi }_{t}u(x)=u(\kappa x+\zeta L),$$

$${\varphi }_{-t}u(x)=u({\kappa }^{-1}x+{\zeta }^{-1}L)\text{ for }-L\le x\le p(\epsilon )\text{ and }\epsilon \le t\le T,$$

$$L_{H_{◃}^1}^2$$

$$L_{(H_{◃}^1{)}^{\ast }}^2$$

$$C_{H_{◃}^1}^k=\left\{\xi \in L_{H_{◃}^1}^2|{\varphi }_{-(\cdot )}\xi (\cdot )\in C^k\left([\epsilon ,T];H_{◃}^1(\epsilon )\right)\right\}\text{ for }k\in \{0,1,\dots \}.$$

$$\dot{\xi }(t):={\varphi }_t\left(\frac{d}{dt}\left({\varphi }_{-t}\xi (t)\right)\right)\in C_{H_{◃}^1}^0,$$

$$W\left(H_{◃}^1,(H_{◃}^1{)}^{\ast }\right)=\left\{u\in L_{H_{◃}^1}^2|\dot{u}\in L_{(H_{◃}^1{)}^{\ast }}^2\right\}.$$

$$T_{◃}:$$

$$(x,t)\mapsto$$

$$-{\partial }_t{\Phi }_{◃}(x,t)-{\partial }_{xx}{\Phi }_{◃}(x,t)$$

$${\partial }_x{\Phi }_{◃}(-L,t)$$

$${\Phi }_{◃}(p(t),t)$$

$${\Phi }_{◃}(x,T)$$

$$-{\partial }_t{\Phi }_{▹}(x,t)-{\partial }_{xx}{\Phi }_{▹}(x,t)$$

$${\partial }_x{\Phi }_{▹}(L,t)$$

$${\Phi }_{▹}(p(t),t)$$

$${\Phi }_{▹}(x,T)$$

Full text

1 Introduction

In this paper we use techniques developed in the field of data assimilation to predict the dynamics of a nonlinear parabolic free boundary price formation model proposed by Lasry & Lions in [16]. The Lasry-Lions (LL) model describes the price evolution of a single good traded between a large group of buyers and a large group of vendors. The price enters as a free boundary, at which trading takes place. After the realisation of a transaction, buyers and vendors immediately sell or rebuy the good at a shifted price. The shift in the price is due to the previously paid constant transaction costs. The situation detailed above can be described by the following nonlinear parabolic partial differential equation

[TABLE]

The positive part $f^{+}=\max(f,0)$ of the function $f=f(x,t)$ corresponds to the distribution of buyers over the price $x\in\Omega$, the negative part $f^{-}=\min(f,0)$ to the is the vendor distribution over the price. The free boundary $p=p(t)$ corresponds to the price where $f(\cdot,t)=0$, the function $\Lambda$ to the total number of transactions executed at that price. The immediate placement and execution of new bids and orders after the trading event are modelled by the Delta Diracs at the shifted prices $p(t)+a$ and $p(t)-a$, where $a\in\mathbb{R}^{+}$ denotes the transaction costs. Random changes in the buyer and vendor distribution are included by a Laplacian with constant diffusivity $\sigma\in\mathbb{R}^{+}$. We assume that the initial distribution $f_{0}$ satisfies:

[TABLE]

and set $\frac{\sigma^{2}}{2}=1$. System (1.1) can be posed on the positive real line $\Omega=\mathbb{R^{+}}$ or a bounded interval $\Omega=[0,x_{\max}]$, where $x_{\max}$ denotes the maximum price. We will consider (1.1) on the bounded interval only and impose homogeneous Neumann boundary conditions

[TABLE]

to ensure that the total number of buyers and vendors is constant in time.

For convenience, we assume the initial price $p_{0}$ is normalized to [math] and only consider its relative change. Hence we work on the shifted domain $[-L,L]$, where $L=\frac{x_{\max}}{2}$. Altogether we will consider (1.1) with boundary condition (1.3) on $\Omega=[-L,L]$ throughout this manuscript.

The LL model (1.1) was analysed in a series of papers, cf. [11, 19, 8, 4, 5]. Most available results are based on a nonlinear transformation of (1.1), which transforms the problem to the heat equation with nonlinear boundary conditions. This connection provides the main analytical ingredients to study existence and long time behaviour of solutions to (1.1). Lasry and Lions introduced the model on the macroscopic level only, a more detailed microscopic interpretation of the trading process and the respective limit as the number of buyers and vendors tend to infinity was missing. This connection was established by Burger et al., who proved that the original LL model can be derived from a Boltzmann type model as the number of transactions tends to infinity, see [2]. In their approach trading events between buyers and vendors are modelled by “collisions”, which can also be used to describe price dynamics in case of more general trading rules. The connection between the Boltzmann-type price formation model and the LL model (1.1) was further investigated in different asymptotic limits in [3]. The LL and Boltzmann-type price formation models are appealing in many respects, especially in terms of analytical tractability. However the resulting price process is deterministic and does not give any insights into connections between transactions rates, order flows or price volatility. Markowich et al., [18] considered a stochastic extensions of the original LL model. However this extension did not give realistic price dynamics either. Very recently Cont and Müller [10] proposed a stochastic partial differential equation with multiplicative noise, which reproduces statistical properties of real price dynamics.

In this paper we focus on the inverse problem of determining the buyer-vendor distribution given measurements of the price and the transaction rate on a time interval $[0,T]$. This distribution can then be used as an initial value and thus allows us to predict price and transaction rate for $t>T$. More specifically we will investigate the question

Problem I: Given measurements of the price $p(t)$ and the transaction rate $\Lambda(t)$ in some time interval $[0,T]$, is it possible to predict the price for times $t>T$ ?

Our approach is based on an optimal control approach proposed by J-P. Puel, see [20, 21]. It is based on a duality argument, which allows to reconstruct the distribution $f$ at the final time $T$. This is in contrast to standard data assimilation where one tries to recover the initial datum $f_{0}(x)$. We adopt the strategy of Puel et al. and use duality estimates to compute linear functionals of $f(T,x)$. These functionals involve the solution of optimal boundary control problems with PDE constraints. Optimal boundary control problems are well studied in the literature, see e.g. [17, 22, 13]. We will make use of an exact null controllability result for parabolic boundary control problems shown in [7]. Its proof is based on Carleman estimates, a technique commonly used to derive exact controllability results (and also uniqueness for inverse problems), see [23, 14] for details. A possible numerical realisation of Puel’s strategy was presented in [9].

Our contributions to the subject of optimal control for parabolic free boundary problems and data assimilation in price formation models are the following:

• •

We present the first approach to reconstruct the buyer- and vendor distribution from measurements of price and transaction rate (to the author’s knowledge).

• •

We generalise the data assimilation approach of Puel et al., see [20], to free boundary value problems and evolving domains.

• •

We provide stability estimates, which give novel insights into the influence of measurement errors on the price dynamics.

• •

We propose a computational strategy to implement the developed framework numerically.

This paper is organized as follows: The proposed framework is based on several analytic results, which will be presented in Section 2. The data assimilation problem itself is discussed in Section 3. Section 4 is devoted to stability in the presence of measurement errors and we conclude by presenting numerical experiments in Section 5.

2 Preliminary results

In this section we provide analytic tools and results of the forward problem and define the respective adjoint problem, which will be used in the optimal control approach.

The presented results rely on the following assumptions:

1. (A1)

$f_{0}(p_{0})=0,\;f_{0}(x)>0\text{ for }x<p_{0}\text{ and }f_{0}(x)<0\text{ for }x>p_{0}$. 2. (A2)

For every $t\in[0,T]$, there exists a constant $\underline{p}>a$ such that $-L+\underline{p}\leq p(t)\leq L-\underline{p}$.

Assumption (A1) is the necessary compatibility condition for the initial datum $f_{0}$ (which we already stated in (1.2)), while (A2) ensures that the price stays sufficiently far away from the interval boundaries. Note that the restriction on $p(t)$ is not severe in the context of inverse problems: Since we will assume later on that we know measurements of $p(t)$ in some time interval $[0,T]$, we can always chose the domain size $L$ (within realistic bounds) such that the condition $p(t)\in(-L+a,L-a)$ is satisfied. As $p(t)$ is continuous, we also know it will stay in $(-L+a,L-a)$ for some time so that it is safe to predict for $t>T$.

2.1 Nonlinear transformation of the model

We start by discussing the nonlinear transformation which converts (1.1) to a linear heat equation. This connection was exploited in almost all analytic results as well as computational methods. It is based on the fact that the second derivative of the buyer vendor distribution $f$ at the price $p(t)\pm a$ behaves like $\Lambda(t)\delta_{p(t)\pm a}$ Thus, shifting the function by multiples of $\pm a$ and adding them up ’eliminates’ the singularity on the right hand side. More precisely, for $\Omega=\mathbb{R}$, we define

[TABLE]

Then the function $F=F(x,t)$ satisfies the heat equation

[TABLE]

with the transformed initial datum

[TABLE]

Since we consider (1.1) with homogeneous Neumann boundary conditions on the interval $(-L,L)$, the sum in (2.1) is finite. If we assume that the initial price is a multiple of $a$, then the transformed initial condition is given by

[TABLE]

with

[TABLE]

We recall that the solution of the original LL model (1.1) can be computed by

[TABLE]

This back-transformation allows us to deduce the corresponding transformed Neumann boundary conditions

[TABLE]

Remark 2.1**.**

As explained above the non-linear transformation is tailored specifically to (1.1), i.e. the fact that first derivates of solutions to (1.1) at $x=p(t)$ are of the form $\pm\lambda\delta_{p(t)\pm a}$ and thus, when summing up shifted solutions the terms on the right hand side vanish. As soon as the structure of the equation is changed, e.g. by adding new non-linear terms or constraints (that act away from $x=p(t)$), the transformation will no longer work.

2.2 Existence and regularity of the price

In the following we provide additional existence and regularity results for the direct problem. Note that they are not optimal in terms of regularity but sufficient to define all quantities that we shall need in the sequel.

Theorem 2.2** (Existence of $f$, $p(t)$).**

Let $f_{0}\in L^{2}(-L,L)$ and $p_{0}\in(-L+a,L-a)$ satisfy (A1). Then the BVP (1.1) has a global solution conserving the total mass of buyers and vendors iff the zero level set $p=p(t)$ of the solution of (2.2)–(2.4) satisfies $p(t)\in(-L+a,L-a)$ for all $t>0$. Then the free boundary $p(t)$ converges to the stationary price $p_{\infty}\in(-L+a,L-a)$.

Proof.

The proof is mainly based on the definition of the transformation (2.1), see [5] for details. ∎

Note that the stationary price is determined by the initial mass of buyers and vendors as well as the transaction rate $a$. In particular

[TABLE]

where $M^{l}=\int_{-}L^{p_{0}}f_{0}(x)dx$ and $M^{r}=\int_{p_{0}}^{L}f_{0}(x)dx$. The presented analysis of the adjoint and assimilation problem relies on the following regularity result for the price $p=p(t)$.

Lemma 2.3** (Regularity of $p(t)$).**

Let $f_{0}\in L^{2}(-L,L)$ and $p_{0}\in(-L+a,L-a)$ satisfy (A1). Then $p(t)\in C^{1}([\varepsilon,T])$ for $\varepsilon>0$.

Proof.

The results is a direct consequence of the fact that $F(x,t)$ is smooth in space and time for all $t>0$ and of the boundedness of $\Lambda$. Indeed, differentiating the relation $F(p(t),t)=0$ yields

[TABLE]

and therefore

[TABLE]

where the parabolic version of Hopf’s Lemma applied at $x=p(t)$ ensures that $\underline{\Lambda}=\min_{t\in[\varepsilon,T]}\Lambda(t)>0$. ∎

Remark 2.4**.**

The regularity of the price $p$ as well as the buyer-vendor density $f$ at the initial time is crucial to define the transformation between the time-dependent domains $[-L,p(t)]$ and $[L,p(t)]$ and the reference domain $[0,1]$ (see Subsection 2.3). It is also important for the exact controllability results of Theorem 3.3. Therefore we will work the temporal domain $[\varepsilon,T]$ instead of $[0,T]$ for some fixed $\varepsilon>0$ in the following only.

2.3 Evolving spaces and the transformation to fixed domains

A crucial step in the subsequent analysis is the splitting of the domain $\Omega$ into the part left and right of the price $p(t)$ (illustrated in Figure 1). We introduce the domains

[TABLE]

as well as

[TABLE]

Following [1], we define evolving Bochner spaces on these domains. We present the construction for the left domain $\Omega_{\triangleleft}=[-L,p(t)]$ only, since the argument for the right domain is analogous. First denote by $H_{\triangleleft}^{1}(t):=H^{1}((-L,p(t))$ the evolving Hilbert space. Next we define the map $\phi_{t}:H_{\triangleleft}^{1}(\varepsilon)\to H_{\triangleleft}^{1}(t)$ by

[TABLE]

with $\kappa=\frac{p(\varepsilon)+L}{p(t)+L}$ and $\zeta=\frac{p(\varepsilon)-p(t)}{p(t)+L}$ for all $-L\leq x\leq p(t)$ and $\varepsilon\leq t\leq T$. The function $\phi_{t}$ is obviously continuous and reduces to the identity at $t=\varepsilon$. It is also a homeomorphism as its inverse

[TABLE]

is continuous as well. This allows us to introduce the evolving Bochner spaces (as in [1, Definition 2.7])

[TABLE]

and, following again [1], make the identification of $u(t)=(\bar{u}(t),t)$ with $\bar{u}(t)$ for $u\in L_{H_{\triangleleft}^{1}}^{2}$ (and likewise in $L_{(H_{\triangleleft}^{1})^{*}}^{2}$). The space of continuously differentiable functions on evolving Bochner spaces is given by

[TABLE]

Thus we can, as in [1, Definition 2.20], to give a notion of time (material) derivative as

[TABLE]

for any $\xi\in C_{H_{\triangleleft}^{1}}^{1}$. Then we can finally define the space used for the notion of weak solutions, namely

[TABLE]

The definitions of the respective quantities $L_{H_{\triangleright}^{1}}^{2}$, $L_{(H_{\triangleright}^{1})^{*}}^{2}$, $C_{H_{\triangleright}^{1}}^{k}$, and $W\left(H_{\triangleright}^{1},(H_{\triangleright}^{1})^{*}\right)$ are analogous.

While the previous definitions allow us to directly work in a noncylindrical domain, it is sometimes also useful consider the transformation to the fixed domain $Q=[0,1]\times[\varepsilon,T]$. Hence we introduce transformations which map $Q_{\triangleleft}$ and $Q_{\triangleright}$ to $Q$:

[TABLE]

Note that due to assumption (A1), $T_{\triangleleft}$ and $T_{\triangleright}$ are well-defined and that $T_{\triangleright}$ actually flips the domain, i.e. it swaps left and right boundary points.

2.4 Adjoint equations

The next ingredient will be two adjoint equations, posed on the domains $Q_{\triangleleft}$ and $Q_{\triangleright}$, respectively.

Definition 2.5** (Adjoint equations).**

For any $\varepsilon>0$, $\psi_{\triangleleft}\in L^{2}(-L,p(T))$, $\psi_{\triangleright}\in L^{2}(p(T),L)$, $u_{\triangleleft},\,u_{\triangleright}\in L^{2}(\varepsilon,T)$ and $T>0$, we introduce the backward in time adjoint equations

[TABLE]

and

[TABLE]

Applying the existence theory of, e.g. [1], for equations on evolving domains, we obtain the following theorem.

Theorem 2.6**.**

Let $p\in C^{1}([\varepsilon,T])$ be given. Then, for every $\Psi_{\triangleleft}\in L^{2}(\Omega_{\triangleleft})$, $u_{\triangleleft}\in L^{2}(\varepsilon,T)$ and every $\Psi_{\triangleright}\in L^{2}(\Omega_{\triangleright})$, $u_{\triangleright}\in L^{2}(\varepsilon,T)$ there exist unique solutions $\Phi_{\triangleleft}$ and $\Phi_{\triangleright}$ to (2.11) and (2.12), respectively. Furthermore, we have

[TABLE]

With the help of the transformations $T_{\triangleleft}$ and $T_{\triangleright}$, equation (2.11) and (2.12) can be transformed into a generic problem of the form

[TABLE]

For (2.11) we define $(y,t)=T_{\triangleleft}(x,t)$ and compute

[TABLE]

while for (2.12) and $(y,t)=T_{\triangleright}$ we obtain

[TABLE]

Note that in view of Lemma 2.3 and Assumption (A1), the coefficients $a$ and $b$ are (in both cases) continuous and uniformly bounded by

[TABLE]

as there may be points with $p^{\prime}(t)=0$. Thus, standard existence and regularity results for linear diffusion–convection equations on fixed domains, such as [15, Theorem 5.2], can be used to ensure the solvability of (2.15).

3 Data assimilation problem

We now turn to the main part of this paper - the inverse or data assimilation problem 1. In classic data assimilation approaches one would use the measurements of $p=p(t)$ and $\Lambda=\Lambda(t)$ on $[0,T]$ to reconstruct the initial datum $f_{0}(x)$ of (1.1). Here we follow an alternative approach proposed by Puel et al., see [20, 21], and estimate the buyer-vendor distribution at the final time, that is $f(x,T)$ instead. This requires the solution of additional optimal control problems, which are, however, well posed if an appropriate regularisation (penalty) is added.

To use Puel’s strategy in our setting, we will estimate the densities of buyers and of vendors separately (that is on the right and left of the free boundary). The reconstruction is based on the following two duality estimates:

Theorem 3.1**.**

Let $f_{0}\in L^{2}(\Omega)$ satisfying assumption (A1) and let $f\in L^{2}(0,T;H^{1}(\Omega))$ be a solution to (1.1) with corresponding price, $p\in C^{1}([\varepsilon,T])$ satisfying (A2). Furthermore, let

[TABLE]

satisfy (2.11) and (2.12), respectively. Then, the following duality identities

[TABLE]

hold for arbitrary functions $u_{\triangleleft},\,u_{\triangleright}\in L^{2}(0,T)$ and every $\varepsilon>0$.

Proof.

We prove the first estimate only, since the argument for (3.2b) is the same. We have

[TABLE]

where we have used the boundary condition (1.3), $f(p(t),t)=0$ and the definition of $\Lambda$. ∎

Now we will use (3.2a)–(3.2b) to determine $f(x,T)$. Since the choice of $\Psi_{\triangleleft}$ and $\Psi_{\triangleright}$ in (3.2a) and (3.2b) is arbitrary and the last term on the right hand side contains only known (i.e. computed or measured) quantities, we could obtain a linear functional of $f(x,T)$. The only unknowns are the first terms on the respective right hand sides. But since we are free to choose arbitrary boundary data $u_{\triangleleft}$ and $u_{\triangleright}$, this leads to the null–controllability problems for (2.11)–(2.12). Indeed, if we can chose $u_{\triangleleft}$ and $u_{\triangleright}$ such that $\Phi_{\triangleleft}(x,\varepsilon)=0$ and $\Phi_{\triangleright}(x,\varepsilon)=0$, the unknown terms in both orthogonality relations drop out and we can reconstruct $f(x,T)$.

3.1 Optimal control problem

To conduct the strategy outlined above, we have to solve the optimal control problems

[TABLE]

Since the structure of both problems is the same, we will only discuss the first one. To increase readability, we will drop the subscript $\triangleleft$ and write $u$, $\phi$, $\ldots$ instead of $u_{\triangleleft}$, $\phi_{\triangleleft}$ from now on. The next result states that the optimal control problem is indeed exactly null-controllable in the sense of the following definition.

Definition 3.2**.**

We say that problem (3.3) is exactly null-controllable, if for every initial datum $\Psi\in L^{2}(\Omega_{\triangleleft})$ to (2.11) there exists $\bar{u}\in L^{2}(\varepsilon,T)$ such that the solution $\Phi$ to (3.3) with control $u=\bar{u}$ satisfies $\Phi(x,\varepsilon)=0.$

The following exact boundary controllability result is based on [7, Theorem 2.3], slightly extended and adapted to our situation. The theorem reads as follows.

Theorem 3.3** (Exact null-controllability).**

For every $\Psi\in L^{2}(\Omega_{\triangleleft})$, there exists at least one control $u\in L^{2}(\varepsilon,T)$ such that the solutions $\Phi$ of (2.11) satisfies $\Phi(x,\varepsilon)=0$ on $\Omega_{\triangleleft}$. Furthermore, there exists a constant $C$ which depends on $p(t)$, $L$ and $T$ such that

[TABLE]

holds with $\bar{u}$ being the control of minimum $L^{2}$–Norm.

Proof.

The regularity of the price $p$ allows us to transform the problem to a fixed domain using $T_{\triangleleft}$ defined in (2.10). Hence we only consider equations of type (2.15). First we observe that for any positive $\delta>0$, any solution $\Phi$ to (2.15) with initial datum $\Psi$ is, by standard parabolic regularity [12, Chapter 7.1], in $L^{2}(\varepsilon+\delta,T);H^{1}(0,1))$ with the estimate

[TABLE]

Thus, we can assume that already $\Psi\in H^{1}(0,1)$ holds. Since by lemma 2.3, $p\in C^{1}([\varepsilon,T])$ (and thus the coefficients $a$ and $b$ in (2.15) are continuous) we can apply [7, Theorem 2.3] to conclude the requested boundary controllability. The continuity estimate (3.5) then follows by combining (3.6), the respective estimate from [7, Theorem 2.1] for the distributed control problem and a standard trace inequality.

∎

In order to be able to numerically solve the optimal control problem, we introduce the following regularized version

[TABLE]

Standard arguments guarantee the existence of a unique minimizer, see e.g. [22, Section 3.5]. Calculating the derivatives of the corresponding Lagrange functional

[TABLE]

we obtain the first order optimality system

[TABLE]

where $\Phi$ satisfies the adjoint equation (2.11) and the coupling

[TABLE]

The following results examine the convergence of $u$ as $\alpha\to 0$. The proofs are using the same techniques as in [21], yet adapted to our boundary control problem.

Theorem 3.4**.**

For every $\alpha>0$, denote by $(u_{\alpha},\phi_{\alpha})$ the corresponding solution to (3.7). Then we have

[TABLE]

where $\bar{u}$ is the solution to the optimal control problem (3.3) having minimal $L^{2}$-norm and $\Phi_{\alpha}$ and $\bar{\Phi}$ are the solutions to (2.11) with boundary data $u_{\alpha}$ and $\bar{u}$, respectively.

Proof.

By Theorem 3.3, we know that there exists at least one function solving the exact null controllability problem. Thus, the set of all these controls in $L^{2}(\varepsilon,T)$ is nonempty. As it is also convex and closed, there exists a unique $\bar{u}$ having minimal $L^{2}$-norm. Since $u_{\alpha}$ minimizes the functional (3.7) among all function in $L^{2}(\varepsilon,T)$ we have

[TABLE]

which implies the (uniform in $\alpha$) bound

[TABLE]

Thus, we can extract a subsequence, again labelled $u_{\alpha}$ that converges weakly to some $\tilde{u}$ in $L^{2}(\varepsilon,T)$. Using the weak formulation of (2.11) and an Aubin-Lions argument, we see that this is sufficient to obtain the convergence

[TABLE]

and (3.13) implies $\tilde{\Phi}(\varepsilon,x)=0.$ Thus, arguing as in the proof of [21, Theorem 2.12], we can use the fact that $\bar{u}$ has minimal norm as well as the lower semi-continuity of the norm w.r.t weak convergence to obtain that $\tilde{u}=\bar{u}$. This argument also implies norm convergence and the uniqueness of the limit then finally yields

[TABLE]

This also implies $\tilde{\Phi}=\Phi$ which completes the proof. ∎

Remark 3.5**.**

Understanding the optimal control problem (3.3) (or (3.4)) as Tikhonov regularisation, one could ask for convergence rates of $u_{\alpha}$ to $\bar{u}$ as $\alpha\to\infty$. Indeed, such rates could be expected under appropriate source conditions on $\bar{u}$. The interesting point now is to understand the influence of $p(t)$ in the definition of the forward operator in the characterisation of such conditions and also how perturbation in $p$ would influence them. We leave this question for future research.

4 Stability in the presence of measurement errors

Assume we have measurements of two different prices $p_{1}$ and $p_{2}$ as well as two different transaction rates $\Lambda_{1}(t)$ and $\Lambda_{2}(t)$. Can we control the difference in the reconstructions $f_{1}(x,T)$ and $f_{2}(x,T)$ as well as the future predicted prices $p_{1}(t)$ and $p_{2}(t)$ for $t>T$ in terms of these differences? In this section we will give a positive answer to this question based on the following strategy

Estimate the error in the optimal controls $u_{1}$ and $u_{2}$ in terms of the error in $p_{1}$ and $p_{2}$ (Lemma 4.2). 2. 2.

Estimate the error in the respective reconstructions $f_{1}(x,T)$ and $f_{2}(x,T)$ in terms of errors in price and transaction rate (Lemma 4.3). 3. 3.

Use these results to predict errors in the future price (Lemma 4.7).

Note however that for the last point we need to make additional regularity assumptions on the reconstructed final data that do not directly follow from our analysis (see Remark 4.5 for details). We start by assuming

1. (A3)

$\|p_{1}-p_{2}\|_{C^{1}([\varepsilon,T])}\leq\delta_{p}\quad\text{and}\quad\|\Lambda_{1}-\Lambda_{2}\|_{L^{2}((\varepsilon,T))}\leq\delta_{\lambda}.$

W.l.o.g. we only consider the optimality system related to (3.3), i.e. the left part $\Omega_{\triangleleft}=[-L,p(t)]$ and again drop the subscript $\triangleleft$. Moreover, we transform all equations to the unit interval $[0,1]$, so that the optimality system reads as

[TABLE]

with $a(t)$ and $b(t)$ as defined in (2.16). Note that the transformed primal and dual equations are still adjoint to one another, yet now with respect to the scalar product

[TABLE]

Lemma 4.1**.**

Let $\Phi$ and $G$ be the solutions to (4.1a) and (4.1e), respectively. Then we have

[TABLE]

with $C_{1}=C_{1}(\alpha,\underline{p},L,T)$ and $C_{2}=C_{2}(\underline{p},L,T)$.

Proof.

These are standard estimates obtained choosing $\Phi$ and $G$ as test functions in the weak formulation of (4.1a) and (4.1e), respectively. For the first estimate, we additionally used the $L^{2}$-bound (3.14) on the boundary control, which introduced the $\alpha$-dependence in $C_{1}$. ∎

Now we are able to prove stability of the optimal control problem in terms of measurement errors in the price.

Lemma 4.2** (Stability of $u$).**

Consider two different prices $p_{1}(t)$ and $p_{2}(t)$ such that $p_{1}(\varepsilon)=p_{2}(\varepsilon)$ and $\|p_{1}-p_{2}\|_{C^{1}([\varepsilon,T])}\leq\delta_{p}$. Denote by $\Phi_{1}$ and $\Phi_{2}$ and $G_{1}$ and $G_{2}$ the solutions to (4.1a) and (4.1e) with $p=p_{1}$ and $p=p_{2}$, respectively. Then the following stability estimate for the controls $u_{1}$ and $u_{2}$ holds:

[TABLE]

Proof.

For each $p_{i}$ (and corresponding $a_{i}$, $b_{i}$), we denote by $G_{i}$, $\Phi_{i}$ and $u_{i}$ the corresponding solutions to the optimality system (4.1a)–(4.1i) and furthermore

[TABLE]

Then, $\bar{\Phi}$ and $\bar{G}$ satisfy, in the weak sense, the equations

[TABLE]

Note that the following calculations are formal since for now we only know existence of weak solutions and therefore some of the integrals are not defined. In the end we arrive, however, at an estimate which is again well defined and could can be obtained rigorously by directly working with weak solutions. We chose this way of presentation as we believe it to be easier to follow. Thus (formally) taking equation (4.3a) and testing it with $\bar{G}$ (with respect to the scalar product (4.2)) yields

[TABLE]

Integrating by parts on the left hand side, using (4.3e) and the boundary conditions results in

[TABLE]

A final integration by parts to remove the second derivatives on the right hand side gives

[TABLE]

Using the estimates of Lemma 4.1, the boundedness of $u$ in $L^{2}$ (see (3.14)) and Cauchy’s inequality applied to the last term on the right hand side, we have

[TABLE]

where we also used the lower bounds (2.18) on $a$ and Assumption (A3) to estimate the expression $(p(t)+L)$ from below by $\underline{p}$ and above by $L-\underline{p}$. Using again (2.18) yields

[TABLE]

Combining this with the previous estimate yields the assertion. ∎

For the second step of our strategy, we return to the orthogonality relation (3.2a) which, transformed to $[0,1]$, reads as

[TABLE]

In the presence of errors in $p$ and $\Lambda$ we obtain two different relations and the following stability result. Note that the above results on the adjoint equations imply solvability for $\Phi$ with continuous dependence on the initial value for any $\Psi\in L^{2}([0,1])$. Hence, the duality relation uniquely defines $f(\cdot,T)\in L^{2}([0,1])$ when given $f(\cdot,\epsilon)\in L^{2}([0,1])$. There is further stable dependence of $f(\cdot,T)$ on the errors in the price and transaction rates, which we make precise by the following result:

Lemma 4.3** (Stability of $f(x,T)$).**

Let $p_{1},\,p_{2}$ and $\Lambda_{1},\,\Lambda_{2}$ be given functions which satisfy Assumption (A3) and denote by $f_{1}(x,T)$ and $f_{2}(x,T)$ the corresponding reconstructed prices calculated using (4.4). Then we have

[TABLE]

Proof.

Subtracting (4.4) for $(p_{1},\lambda_{1},u_{1})$ and $(p_{2},\lambda_{2},u_{2})$ yields

[TABLE]

We estimate each term of the right hand side separately

[TABLE]

where we used Lemmata 4.1 and 4.2. Next we have

[TABLE]

using that for positive times $t\geq\varepsilon$ (and away from the boundary) $\phi_{1}$ is Lipschitz continuous. Next we have

[TABLE]

and finally

[TABLE]

Combining all estimates and taking the supremum over all $\Psi\in L^{2}((0,1))$ with $\|\Psi\|_{L^{2}((0,1))}=1$, we finally obtain

[TABLE]

Taking $C_{6}=\max(C_{13},C_{11})$ yields the assertion. ∎

Remark 4.4**.**

The estimates of Lemma 4.2 and 4.3 show that, for $\alpha>0$, the reconstruction of the unknown buyer vendor distribution $f(x,T)$ is actually a well-posed problem, at least for sufficiently smooth perturbations of $p$. This is due to the fact that we are solving a regularized optimization problem. The price to pay is that the term involving $f(x,\varepsilon)$ in (4.4) does not vanish. However, since $f(x,\varepsilon)$ is fixed, is does not appear in our stability estimates.

For the next result, we choose perturbed prices $p_{1}$ and $p_{2}$ such that $|p_{1}(T)-p_{2}(T)|<2a$ and assume w.l.o.g. that $p_{1}(T)\leq p_{2}(T)$ and make the following additional assumptions:

1. (A4)

$\delta_{p}<a$, 2. (A5)

$f_{1}(x,T),\,f_{2}(x,T)\in L^{2}(-L,L)\cap H^{4}(I)$ with $I\subset(p_{2}(T)-a,p_{1}(T)+a)$ 3. (A6)

$\|f_{1}(x,T)-f_{2}(x,T)\|_{H^{4}(I)}\leq C_{6}(\delta_{p}+\delta_{\lambda})$

Remark 4.5**.**

We mention that indeed it is natural to assume strong regularity of $f$ in a neighbourhood of $p(T)$ for $T>0$, since it locally arises as the solution of a heat equation. On the other hand we need to expect some singularities around $p(T)-a$ and $p(T)+a$ due to the singular source terms. Thus (A5) seems completely natural for forwards solutions of the price formation model. Moreover, it can also be verified that $f(\cdot,T)$ reconstructed via (4.4) has local $H^{4}$-regularity, which follows from using $\Psi$ supported in $I$ and an analysis of the solution of the parabolic equation for $\Phi$, which can be estimated in terms of the $H^{-4}$ norm of the initial value.

In the following we analyse the forward propagation for $t>T$ in a small time interval. We denote the new initial value by $f_{i,0}:=f_{i}(\cdot,T)$. First note that using the same localisation strategy as in [19] (i.e. multiplying the solution to (1.1) with a smooth cut-off function that has support inside the interval $I$), implies

[TABLE]

with $\gamma>0$ to be fixed later on and where $f_{i}$ is the solution to (1.1) with the reconstructed initial datum $f_{i}(x,T)$ that additionally satisfies (A4)-(A6). Furthermore, $I_{2}$ is an interval that is compactly contained in $I$. This allows us to derive the following estimates on terms of the form $(\partial_{x}K_{N})\ast f$, where we denote by $K_{N}(x,t)$ the heat kernel with Neumann boundary conditions on $[-L,L]$, see e.g. [6, Section 6.4], and furthermore use the notation

[TABLE]

Lemma 4.6**.**

For given $T>0$, initial values $f_{0},\,f_{1,0}$ and $f_{2,0}$ at time $T$ satisfying (A4), (A5) we have for $t\in[T,T+\gamma]$ with $\gamma$ sufficiently small

[TABLE]

For two continuous functions $p_{1}(t)$ and $p_{2}(t)$, we have

[TABLE]

Proof.

First note that $(\partial_{x}K_{N}^{T})\ast f$ is the solution to the heat equation with homogeneous Neumann boundary condition, zero right hand side and initial datum $f$. Then, the first estimate is a direct consequence of (4.6) applied to such an solution with initial datum $f_{1,0}-f_{2,0}$. The second one follows from the fact that, as for $t$ sufficiently small, $p_{1}(t)$ and $p_{2}(t)$ are in $I_{2}$ and thus, using again (4.6), the derivative of a solution to the heat equation that appears on the left hand side is Lipschitz continuous. ∎

We are now in a position to state the stability result for future prices.

Lemma 4.7**.**

Let assumptions (A3)–(A6) be satisfied and denote by $f_{1}$ and $f_{2}$ the solution to (1.1) on the time interval $[T,T+\gamma]$ with initial data $f_{1}(x,T)$ and $f_{2}(x,T)$, reconstructed from measurements $p_{1},\Lambda_{1}$ and $p_{2},\Lambda_{2}$ in $[0,T]$, respectively. Then there exists a constant $\gamma>0$ and we the corresponding prices $p_{1}$ and $p_{2}$ for $t\in(T,T+\gamma)$ satisfy the estimate

[TABLE]

Note that unfortunately, we cannot give a lower bound on the quantity $\gamma$ as it depends in a non-linear and non-local fashion on the initial datum via the solution of the equation. However, the proof below shows that as the transaction rate $\Lambda$ increases, also $\gamma$ becomes larger which agrees with the modelling.

Proof.

Due to assumption (A5) we can invoke [19][Lemma 2.5] to show that for $\gamma$ sufficiently small (depending on $f_{i}(x,T)$, $i=1,2$) the corresponding transaction rates $\Lambda_{1},\,\Lambda_{2}$ are strictly positive on $[T,T+\gamma]$. Furthermore, (A5) implies that $p_{1}(t),\;p_{2}(t)$ are in $C^{1}([T,T+\gamma])$. Now Duhamel’s formula allows us to express the solutions $f_{i}(x,t)$ to (1.1) as

[TABLE]

Taking the space derivative and evaluating at $x=p_{i}(t)$ we obtain

[TABLE]

Subtracting (4.8) for $i=1,2$ and using the linearity of the convolution, we obtain

[TABLE]

with $\theta_{i}(t,\tau)=[\partial_{x}K_{N}^{T}(p_{i}(t)-p_{i}(\tau)+a,t-\tau)-\partial_{x}K_{N}(p_{i}(t)-p_{i}(\tau)-a,t-\tau)].$ As the $p_{i}(t)$ are continuous, choosing $\gamma$ sufficiently small guarantees that the derivatives of $K_{N}^{T}$ appearing in the definition of $\theta_{i}$ are always evaluated away from their singularity, in particular they are bounded and locally Lipschitz-continuous, which implies with the local Lipschitz constant $\lambda$

[TABLE]

Taking the absolute value on both sides of (4.9) and using Lemma 4.6 implies

[TABLE]

so that Gronwall’s lemma implies, together with (A3) and (A6), yields

[TABLE]

Next we exploit the fact that $f_{i}(p_{i}(t),t)=0$ by taking the time derivative, which gives

[TABLE]

Subtracting the above equation for $i=1$ and $i=2$ respectively, using the definition of $\Lambda_{i}$ and integrating in time we obtain, for $T\leq t\leq T+\gamma$

[TABLE]

Denoting by $\underline{\Lambda}=\inf_{T\leq s\leq T+\gamma}\Lambda_{1}(s)\Lambda_{2}(s)$ and using (A3) this yields

[TABLE]

As a consequence of (4.6), $\partial_{t}f_{i}(\cdot,t)$ is bounded and Lipschitz continuous. Thus using (4.10), (A3) and once more (4.6) applied to $\partial_{t}f_{1}(p_{1}(t),t)-\partial_{t}f_{2}(p_{1}(t),t)$ (and together with (A6)) finally yields the assertion. ∎

5 Numerical Simulation

We conclude by illustrating the proposed methodologies and confirming the obtained analytic results with various computational experiments. All simulations are performed on the domain $[-L,L]$, which is split into $N$ intervals of length $h$. The discrete grid points are denoted by $x_{i}=ih$. We compute solutions at discrete times $t^{k}=k\Delta t$, where $\Delta t$ is the discrete time step. However we will omit all full time-discrete expressions in the following, to enhance readability.

The reconstruction of the buyer-vendor distribution is based on piecewise linear basis functions. Let $V_{h}$ denote the space of piecewise linear basis functions $\phi_{j}$, which satisfy $\phi_{j}(x_{i})=\delta_{ij}$. We wish to reconstruct $\hat{f}\in V_{h}$, which is given by $\hat{f}(x,T)=\sum_{j=1}^{J}\hat{f}_{j}\phi_{j}(x)$ using the duality estimates (3.2).

Data generation: We solve the transformed LL model (2.2) for a given initial buyer-vendor distribution $f_{0}$. In doing so we transform the initial distribution $f_{0}$ via (2.1), and compute the solution to the heat equation (2.2) using an implicit in time discretization. The returned discrete price $p^{(k)}=p(t^{k})$ corresponds to the zero levelset of the buyer-vendor distribution $F^{(k)}=F(t^{k})$ (computed via linear interpolation). Note that we use a finer spatial and temporal discretization to generate the data than in the subsequent reconstruction.

Steepest descent: We solve (3.7) and the corresponding problem on $\Omega_{\triangleright}$ and $\Omega_{\triangleleft}$ using steepest descent. In doing so, we compute the variational derivatives of (3.8) and obtain the first order optimality system (3.9) as well as the updates for the controls $u_{1}$ and $u_{2}$. The detailed steps are outlined in the While-Loop of Algorithm 5.1. Here the parameter $\beta^{(l)}>0$ is the step size of the steepest descent update. We use the Armijo-Goldstein condition to adjust $\beta^{(l)}$ in the search direction $\mathbf{p}$. We recall that the Armijo-Goldstein condition for a general functional $\mathcal{J}$ is given by

[TABLE]

where $\gamma\in(0,1]$. The starting value is set to $\beta^{(0)}=0.25$, which is then reduced (up to a maximum of four times) by a factor of $\frac{1}{2}$ until condition (5.1) is satisfied. Note that we transform the computational domains $\Omega_{\triangleleft}$ and $\Omega_{\triangleright}$ to $[0,1]$ as discussed in Section 2.3 in all simulations. We solve the forward as well as the adjoint equations using an implicit in time discretization and piecewise linear basis functions in space.

Identifiabilty for different initial conditions.

In the first experiment we set $L=0.5$ and the final time to $T=0.25$. We split the spatial domain $[-0.5,0.5]$ into $200$ elements and the time interval $[0,0.25]$ into $125$ time steps. The initial datum is set to

[TABLE]

We approximate the final buyer-vendor distribution using $J=50$ basis functions. Furthermore we choose the following parameters

[TABLE]

Figure 2 shows the reconstructed and computed function $F$ (the latter computed by solving the heat equation (2.2) with the transformed initial datum $F_{0}$). We observe a good agreement, with small artefacts at the boundary and the buyer-vendor interface. The corresponding controls are shown in Figure 3.

Next we choose a slightly different initial datum, in particular

[TABLE]

In this case the price is not monotone, see Figure 4. However, the quality of the reconstructions is comparable to the one of the previous example.

Stability of $\hat{f}(x,T)$

Next we are interested in the stability of the reconstruction $\hat{f}(x,T)$ with respect to perturbations in the price. Lemma 4.3 and in particular estimate (4.5) state that the difference in the reconstructions is bounded by the difference in the prices and transaction rates. We consider the following perturbation of the unperturbed price $p_{\delta,0}$:

[TABLE]

Note that the perturbed price is still in $C^{1}$ and that $p_{\delta,0}(0)=p_{\delta,k}(0)$ for all $k=1,\ldots K$. We start with the same the initial datum as in the first example, that is (5.2), and set $\delta=0.01$ and $K=13$. All other parameters are the same as in the first example, except that we reconstruct the final profile using $J=80$ basis functions. Figure 5 illustrates the linear increase of the error in the controls and the reconstruction as the noise level increases. For larger values of $\delta$ this agrees with the theoretical results of Lemma 4.2 and Lemma 4.3, respectively. For small $\delta$ a saturation effect due to the numerical discretization error occurs.

Predicting price dynamics

We conclude with an example where we use the reconstructed buyer-vendor distribution $f(x,T)$ to estimate future price dynamics and illustrate the influence of noise in those. We consider an initial datum of the form

[TABLE]

on the domain $\Omega=[-0.5,0.5]$, that is $L=0.5$. The computational domain is split into $100$ intervals, the time horizon $[0,0.5]$ into $100$ time steps. Since the number of buyers equals the number of vendors the exact price is constant and equal to zero. The computed price curve, that is the dark blue line in Figure 6, converges quickly to the true equilibrium price (being slightly off at the beginning due to numerical errors caused by the nonlinear transformation of the initial datum).

We then use the original price curve on the time interval $[0,0.5]$, denoted by $p_{\delta,0}$, as well as the perturbed prices

[TABLE]

to reconstruct the solution at time $t=0.5$, that is $f(x,T=0.5)$. Note that the sinusoidal perturbations ensure that the final distribution of buyers and vendors equals the initial distribution. These reconstructed buyer-vendor profiles are then used to compute the future price dynamics. We observe a jump in the price in all prices at time $t=0.5$, which is caused by the numerical error made in the reconstruction. However, the price then converges quickly to an equilibrium value, which is close to the true price (dark blue line).

The parameters used in the reconstruction Algorithm 5.1 are set to

[TABLE]

The final buyer vendor distribution $f(x,T=0.5)$ was computed using $J=80$ basis functions.

6 Summary and Outlook

We studied a data assimilation problem for a parabolic nonlinear free boundary problem. This partial differential equation describes the evolution of the price, that is the free boundary, in a large economic market. We developed an analytical and computational framework for the corresponding data assimilation problem, which is based on a previous work by Puel et al., see [20]. The free boundary splits the original problem into two parts, each of them defining a separate optimal control problem. We discussed analytic properties of the respective problems and derived stability estimates for the controls and reconstructed unknown buyer-vendor distribution in the presence of noise. Finally we confirmed and illustrated our results with computational experiments.

We believe that the developed framework provides the basis for more general data assimilation problems in price formation. In [2] Burger et al. considered a Boltzmann type price formation model, which allows for more complex trading mechanisms. This problem is a system of nonlocal reaction-diffusion equations on the whole domain, where multiple prices (even with continuous distribution) and transaction rates can appear. Analogous questions can be asked for this problem if only the expectation of the price is to be predicted, but the problem could also be extended to a stochastic distribution of the price.

Acknowledgements

MB acknowledges support by ERC via Grant EU FP 7 - ERC Consolidator Grant 615216 LifeInverse. MTW acknowledges financial support from the Austrian Academy of Sciences ÖAW via the New Frontiers Grant NST-001 and the EPSRC via the First Grant EP/P01240X/1.