Construction of a Class of Copula Using the Finite Difference Method

The definition of a copula function and the study of its properties are at the same time not obvious tasks, as there is no general method for constructing them. In this paper, we present a method that allows us to obtain a class of copula as a solution to a boundary value problem. For this, we use the finite difference method which is a common technique for finding approximate solutions of partial differential equations which consists in solving a system of relations (numerical scheme) linking the values of the unknown functions at certain points sufficiently close to each other.


Introduction
Copulas thus appear to be a natural tool for constructing multivariate distributions via Sklar's theorem when the marginal distributions are sufficiently regular. Sklar's theorem thus provides a canonical representation of a multivariate distribution, via the data of marginal distributions and dependency structure. The definition and properties of copulas are illustrated in [1]. In this paper, we will limit ourselves to two variables (dimension 2) for the sake of clarity and conciseness.
The objective of this paper is to construct a family of copulas C, solution of an elliptic partial differential equation (PDE), using the finite difference method. Thus, this new class of copulas solution of the boundary value problem is given by the following problem P ð Þ: where k ≥ 0 a real, C and f are the two copulas, and Δ is the Laplacian operator. The boundary value problem described above can be thought of as a nonhomogenous Dirichlet problem. The solution to this problem is not analytically known in general. An approximation is then made to reduce the problem to a finite number of unknowns (discretization process).
For this, we use the finite difference method which is one of the oldest methods of numerical simulation which is still used for some applications. This method appears to be the simplest to implement because it proceeds in two steps: on the one hand the discretization by finite differences of the differentiation operators and on the other hand the convergence of the numerical scheme thus obtained when the distance between the points decreases. Finite differences are particularly suitable for rectangular or parallelepipedic domains.
We therefore introduce a mesh of steps h in the direction u and v. The nodes of the mesh are the points P ij = ðu i , v j Þ where the solution is approached. The following can be noted: (1) u i = ih, 0 ≤ i ≤ N + 1, the vertices of the mesh in the u direction (2) v j = jh, 0 ≤ j ≤ N + 1, the vertices of the mesh in the v direction We are looking for an approximation of the equation at the nodes of the mesh ðP ij , 1 ≤ i ≤ N ; 1 ≤ j ≤ NÞ. The principle of the finite difference method consists in approximating the derivatives of a function by linear combinations of the values of this function at the points of the mesh. The difficulty of this method is that the solution is only defined by points, so the reconstruction at arbitrary locations is not uniquely defined and the boundary conditions tend to be complicated to implement.

Materials and Methods
In the rest of the document, we will note I = ½0, 1 and I 2 = I × I = ½0, 1 × ½0, 1.
Definition 1 (see [1]). A copula is a function C : I 2 ⟶ I with the following properties: where B is the rectangle ½u 1 , u 2 × ½v 1 , v 2 and the expression (3) defines the C-volume of B.
According to Sklar's theorem (see [1], Theorem 2.3.3) with given F, G, and H defined as above, there exists a copula C such that for all , y ∈ ℝ: and conversely, for any copula C, the function H defined with Equation (2) is a joint distribution function with margins FðxÞ and GðyÞ.The aim of this paper is to construct a family of copula C using the concept of finite difference schemes.
Definition 2 (consistency error). We called error of consistency of the numerical scheme A h u h = b h the vector ε h ðuÞ ∈ ℝ N defined by the following: : where π h ðuÞ represents the projection of the exact solution onto the mesh. The scheme is said to be consistent for the norm k:k of ℝ N if lim h⟶0 ε h ðuÞ = 0. If moreover there is a constant α independent of h such that the schema is said to be of order p for the normk:k.
Usually, we use the norms k:k = k:k 1 , k:k 2 or k:k = k:k ∞ .
Definition 3 (stability). We will say that a scheme is stable for the norm k:k ∞ , there exists a constant α > 0 independant of h such that Theorem 4 (Lax's theorem). A scheme is convergent if and only if this scheme is consistent and stable.

Results and Discussions
In this section, we will solve problem ðP Þ using finite difference schemes. We will observe two cases: when k = 0 and when k > 0.
3.1. Solving ðP Þ When k > 0. In this section, we approach problem ðP Þ with two finite difference schemes. The consistency and possible stability of these schemes will be studied.
In the rest of the work, we will use the scheme which converges the fastest (for simulation with MATLAB). Thus, the following first result is obtained using five points of discretization.
ðP Þ can be written as follows: Journal of Function Spaces where α 0 = 4 h 2 , which are Scheme 1or which are Scheme 2.
The points inside the grid must check the following condition: Proof. Each derivative is discretized according to its own direction, so by applying Taylor's formula in the u and v directions, we have the following: By summing up Equations (10) and (11), we get the following: Equation (12) can be approached as follows: By a similar reasoning, we obtain that Let us now write the second-order approximation of ∂Cðu, vÞ/∂u. By performing the subtraction between Equations (10) and (11), we obtain the following: The first-order approximation of ∂Cðu, vÞ/∂u is obtained by writing the following: then By adding Equations (13), (14), and (15), we have the following: Finally, we get By posing α 0 = 4/h 2 , α 1 = 1/h 2 + k/2h, α 2 = 1/h 2 − k/2h, and α 3 = α 4 = 1/h 2 , we get Scheme 1.
Proof. Let us pose For Scheme 1, we can write the following:

Journal of Function Spaces
As C ∈ C 4 ðIÞ, there exist ε i,j ∈ I and η i,j ∈ I such that Analogous formulas are obtained for C i,j+1 and C i,j−1 , and the following is concluded: where k:k ∞ denotes the uniform norm (see [2]) on I of the derivative h 2 of C relative to u. We finally get where α depends only on C and k. As for h small, we have h 2 ≤ h; it can be deduced that Scheme 1 is consistent. ☐ For Scheme 2, we have the following: and we deduce Finally, where β depends only on C and k. it can be deduced that the Scheme 2 is consistent.
We can only conclude that there is a contradiction if ð1/h 2 − k/2hÞ ≥ 0 ⟹ h ≤ 2/k. Scheme 1 verifies the following relation: when h verifies the condition which will be called stability condition: hk ≤ 2.
Equation (29) shows us the stability of Scheme 1.
Remark 8. Both schemes are convergent: Scheme 1 is of order 2 (in h 2 ) and Scheme 2 is of order 1 (in h). Scheme 1 converges faster but has a condition of stability k ≤ 2/h. Scheme 2 is unconditionally stable.

3.2.
Solving ðP Þ When k = 0. When k = 0, we are dealing with Poisson's equation. Let us suppose the following practical case: the temperature distribution Cðu, vÞ at the ðu, vÞ point of a plate whose sides have a temperature imposed on the edge and which receive an external heat input of density f is modelled by a partial differential equation.

Journal of Function Spaces
ðP Þ can be written as follows: where α 0 = 4 h 2 , The points inside the grid must check the following condition: In addition, Equation (30) admits a unique solution.
Proof. When k =0, we get the two-dimensional Poisson equation. Demonstrations obtain in the same way as in the previous section (by setting k = 0). The obtained diagram admits a single solution and converges. ☐

Numerical Simulation Using MATLAB
We make in this part a simulation of the approached solution using the MATLAB environment. We remind you that we are interested in a temperature diffusion on a plate. This temperature distribution is noted Cðu, vÞ at the ðu, vÞ point and generated by an external source noted f ðu, vÞ at the ðu, vÞ point. The source f being a copula, we give an approximation of the solution C of problem P, using different types of copulas (concerning the copula f ). The method used for the simulation is that of Gauss-Seidel. We will use Scheme 1 because it converges faster than Scheme 2.
For these reasons, we will have the following: (1) The numerical solution of the approximate copula C by finite difference (2) The number of iterations of the Gauss-Seidel method 4.1. Approached Solution of the Problem ðP Þ. Using the MATLAB software, we give an approximate value of the copula. We still have to choose the copulas that will be the source functions.
Let us suppose that f coincides with the independent copula defined by f ðu, vÞ = uv such that u, v ∈ I. We obtain the following figures: From Figure 1, we notice that the figure on the left gives us the number of iterations for Scheme 1 to be stable at 51 iterations.
Let us suppose that f is an Archimedean copula (see [1]). Let us choose two Archimedean copulas: the first is the parametric family of copulas known as the Gumbel-Hougaard family and the second is known as the Ali-Mikhail-Haq family.
The expression of Gumbel-Hougaard family is given by the following: and the Ali-Mikhail-Haq family is given by the following: (i) When f is the Gumbel-Hougaard family, we get the following figure: For Figure 2, the number of iterations required for stability is 960 iterations.  4.2. Case Where the Scheme Is Unstable. Recall, according to Proposition 7, Scheme 1 is stable only if hk ≤ 2. In this section, let us take an example where this is not the case. The goal is to observe the behavior of the infinite norm as a function of the number of iterations. We will take k = 70 and N = 15; then, we have hk = 4:3750 > 2; it follows the following figure: From Figure 6, we clearly notice that the infinite norm of the copula is not bounded, i.e., |kCk ∞ ⟶ +∞. The scheme is so unstable in this case.

Conclusion
In this paper, we constructed a new class of copula by the finite difference method. This method gives an approached solution which converges to the exact solution of the boundary value problem. We approached the problem with two finite difference schemes and chose the fastest converging one. We use MATLAB environment to make numerical simulations.

Data Availability
No data were used to support this study because all figures are made by simulation on MATLAB.

Conflicts of Interest
The authors declare that there are no competing interests regarding the publication of this paper.