Dynamics of a Bertrand Duopoly Game of the Greek Oil Market and Application of the d-Backtest Method

Georges Sarafopoulos¹, Evaggelos Drimpetas¹, Kosmas Papadopoulos^1*, Dimitrios Vezeris²

¹Department of Economics, Democritus University of Thrace, Komotini, Greece.
²Department of Electrical & Computer Engineering, Democritus University of Thrace, Xanthi, Greece.

DOI: 10.4236/am.2021.1211070                

Abstract

This work presents the complexity that emerges in a Bertrand duopoly between two companies in the Greek oil market, one of which is semi-public and the other is private. The game uses linear demand functions for differentiated products from the existing literature and asymmetric cost functions that arose after approaches using the published financial reports of the two oil companies (Hellenic Petroleum and Motor Oil). The game is based on the assumption of homogeneous players who are characterized by bounded rationality and follow an adjustment mechanism. The players’ decisions for each time period are expressed by two difference equations. A dynamical analysis of the game’s discrete dynamical system is made by finding the equilibrium positions and studying their stability. Numerical simulations include bifurcation diagrams and strange attractors. Lyapunov numbers’ graphs and sensitivity analysis in initial conditions prove the algebraic results and reveal the complexity and chaotic behavior of the system focusing on the two parameters k₁ and k₂ (speed of adjustment for each player). The d-Backtest method is applied through which an attempt is made to control the chaos that appears outside the stability space in order to return to the locally asymptotically stable Nash equilibrium for the system.

Keywords

Duopoly Game, Discrete Dynamical System, Nash Equilibrium, Stability, Bifurcation, Chaotic Behavior, d-Backtest Method

1. Introduction

The Cournot’s game was modified in 1883 by the French mathematician Joseph Louis Francois Bertrand. He suggested that the players (sellers) actually choose prices rather than the production quantities. The basic premise of the Bertrand model originally is that all players make decisions in a naive way so that in every step, each player assumes the last values were taken by the competitors without an estimation of their future reactions. However, such an assumption is very unlikely under the conditions of the real market, since not all players share naive beliefs. Some of the authors proposed different approaches to firm behavior, considering duopolies under homogeneous expectations and found a variety of complex dynamics in their games, such as the appearance of strange attractors [1] - [11]. Also, models with heterogeneous agents were studied [2] [12] - [23].

The producers in the real market behave without knowing the entire demand function, though it is possible that they have a perfect knowledge of technology, represented by the cost function. Here it is more likely that firms employ some local estimate of the demand. This issue has been previously analyzed [24] - [31]. Bounded rational players (sellers) update their strategies based on discrete-time periods using a local estimate of the marginal profit. With such a local adjustment mechanism, the players are not requested to have a complete knowledge of the demand and cost functions [11] [13] [32] [33].

The present paper is a partial approach to the authors’ main ongoing research objective, which is the emergence of complexity in various oligopoly models as well as its control. In this study, the dynamics of a Bertrand-type duopoly game with differentiated goods where both firms behave with homogeneous expectation strategies with asymmetric cost functions that are resulted from the regression analysis using their published financial reports. In addition to the differentiation of the cost functions for the two players, there is another assumption that the composition of two companies is different and more specifically that the first one is a private company and the second is public by a percentage (semi-public). It is shown that the model gives more complex, chaotic and unpredictable trajectories as a consequence of the change in the bounded rational players’ speed of adjustment (parameter k). The paper is organized as follows: in Section 2, the dynamics of the duopoly game with homogeneous expectations, differentiated products, linear demand and asymmetric cost functions for two players are analyzed. Both players are set as bounded rationality. The existence and the local stability of the equilibrium positions are studied; in Section 3, numerical simulations are used to verify the algebraic results of Section 2, we show the complex dynamics via graphing strange attractors, computing Lyapunov numbers, and sensitive dependence on initial conditions. Finally, in Section 4, the application of the d-Backtest in this duopoly game is presented as an attempt to control the chaos of the discrete dynamical system that appears.

2. The Game

2.1. The Construction of the Game

The game is based on the assumption that two homogeneous players participate in it and more specifically it is considered that they decide their products’ prices in a rational way following an adjustment mechanism (bounded rational players). They produce differentiated goods and offer them on a common market at discrete time periods $t=1,2,\cdots$. At each period t, each company must form an expectation of the rival’s strategy in the next time period in order to determine the corresponding profit-maximizing prices for periodt + 1. With q₁, q₂ are symbolized the production quantities of two players. It is considered that the preferences of consumers represented by the equation:

$U\left(q_1,q_2\right)=\alpha \left(q_1+q_2\right)-\frac{1}{2}\left(q_1^2+q_2^2+2\cdot d\cdot q_1\cdot q_1\right)$ (1)

where α > 0 is a positive parameter which expresses the market size and $d\in \left(-1,1\right)$ is the parameter that reveals the differentiation degree between two players’ products. Explaining the values of the differentiation parameter d, when it takes the maximum positive value and d = 1, then each product is a substitute product for the other, since the products are homogeneous. It is understood that for positive values of the parameter d, the larger the value, the less diversification there is between two products. On the other hand, negative values of this parameter describe that both products are complimentary and when d = –1, then the phenomenon of full competition between the two companies is appeared. The inverse demand functions (as functions of quantities) coming from the maximizing of Equation (1) are given by the following equations:

$p_1\left(q_1,q_2\right)=\alpha -q_1-d\cdot q_2$ and $p_2\left(q_1,q_2\right)=\alpha -q_2-d\cdot q_1$ (2)

Calculating the direct demand functions (as functions of prices p₁ and p₂) it gives the equations:

$q_1\left(p_1,p_2\right)=\frac{\alpha \left(1-d\right)-p_1+d\cdot p_2}{1-d^2}$ and $q_2\left(p_1,p_2\right)=\frac{\alpha \left(1-d\right)-p_2+d\cdot p_1}{1-d^2}$ (3)

In this work, data were collected from the financial reports of two companies. The data related to sales quantities in thousands of metric tons (k∙MT) and the total revenue in millions of euro (m∙? for each quarter of the years from 2011 to 2020. With this data and using regression analysis an attempt is made to approach the form of the cost function of each company. According to this analysis it can be assumed that the cost function of the first player is quadratic and of the second player is linear¹. The cost functions approaches are shown in Figure 1 (Motor Oil Group) and Figure 2 (Hellenic Petroleum Group).

It is supposed that the cost function of the first player (Motor Oil Group) is:

$C_1\left(q_1\right)=c_1\cdot q_1^2+c_2\cdot q_1+c_3$ (4)

with $c_1>0,c_2\ge 0,c_3\ge 0$ and $c_2^2\le 4\cdot c_1\cdot c_3$ (quadratic cost’s conditions) and the cost function of the second player (Hellenic Petroleum Group) is supposed as linear:

Figure 1. Quadratic cost function approaching of the player 1 (motor oil group).

Figure 2. Linear cost function approaching of the player 2 (Hellenic petroleum group).

$C_2\left(q_2\right)=c_4\cdot q_2+c_5$ (5)

where $c_4>0$ and $c_5\ge 0$ (linear cost’s conditions).

With these assumptions the profits of the two firms are given by:

$\begin{array}{l}{\Pi }_1\left(p_1,p_2\right)=p_1\cdot q_1-C_1\left(q_1\right)\\ =\left(p_1-c_2\right)\cdot \frac{\alpha \left(1-d\right)-p_1+d\cdot p_2}{1-d^2}-c_1\cdot {\left[\frac{\alpha \left(1-d\right)-p_1+d\cdot p_2}{1-d^2}\right]}^2-c_3\end{array}$ (6)

and:

${\Pi }_2\left(p_1,p_2\right)=p_2\cdot q_2-C_2\left(q_2\right)=\left(p_2-c_4\right)\cdot \frac{\alpha \left(1-d\right)-p_2+d\cdot p_1}{1-d^2}-c_5$ (7)

Then the marginal profits at the point of the strategy space are given by:

$\begin{array}{l}\frac{\partial {\Pi }_1}{\partial p_1}=\frac{1}{{\left(1-d^2\right)}^2}[\alpha \cdot \left(1-d\right)\cdot \left(1-d^2+2c_1\right)+c_2\left(1-d^2\right)\\ -2\left(1-d^2+c_1\right)\cdot p_1+d\cdot \left(1-d^2+2c_1\right)\cdot p_2]\end{array}$ (8)

and:

$\frac{\partial {\Pi }_2}{\partial p_2}=\frac{\alpha \cdot \left(1-d\right)+c_4-2p_2+d\cdot p_1}{1-d^2}$ (9)

The Consumer Surplus (CS) is calculated, because it is contained in the social welfare that the second player cares about a percentage of it and it is given by the equation:

$CS=U(q_1,q_2)-p_1\cdot q_1-p_2\cdot q_2=\frac{1}{2}{q}_1^2+\frac{1}{2}{q}_2^2+\left(1+d\right)\cdot q_1\cdot q_2$ (10)

The social welfare (W) is given by:

$W\left(q_1,q_2\right)=CS+{\Pi }_1+{\Pi }_2$ (11)

The second player the Hellenic Petroleum Group (semi-public firm) cares about the maximization of a function that contains a percentage combination between the social welfare and his profit function (generalized relative profit function). This utility function V₂ is described by the following equation:

$V_2\left(q_1,q_2\right)=s\cdot W\left(q_1,q_2\right)+\left(1-s\right)\cdot {\Pi }_2\left(q_1,q_2\right)$ (12)

where $s\in \left[0,1\right]$ and the marginal function of V₂ is given by:

$\begin{array}{c}\frac{\partial V_2}{\partial p_2}=s\cdot \left(\frac{\partial CS}{\partial p_2}+\frac{\partial {\Pi }_1}{\partial p_2}\right)+\frac{\partial {\Pi }_2}{\partial p_2}\\ =s\cdot [-\left(1+2c_1\cdot d\right)\cdot \frac{\alpha \cdot \left(1-d\right)-p_1+d\cdot p_2}{{\left(1-d^2\right)}^2}\\ +\left(d^2+d-1\right)\cdot \frac{\alpha \cdot \left(1-d\right)-p_2+d\cdot p_1}{{\left(1-d^2\right)}^2}+\frac{d\cdot \left(p_1-c_2\right)}{1-d^2}]\\ +\frac{\alpha \cdot \left(1-d\right)+c_4-2p_2+d\cdot p_1}{1-d^2}\end{array}$ (13)

The following substitutions are defined to simplify the mathematical calculations:

• $A=\alpha \cdot \left(1-d\right)\cdot \left(1-d^2+2c_1\right)+c_2\cdot \left(1-d^2\right)$

• $B=2\left(1-d^2+c_1\right)$

• $C=d\cdot \left(1-d^2+2c_1\right)$

• $D=\alpha \left(1-d\right)\cdot \left(1-d^2\right)\cdot \left(1-s\right)+\left(1-d^2\right)\cdot \left(c_4-s\cdot c_2\cdot d\right)-2s\cdot \alpha \left(1-d\right)\cdot d\cdot c_1$

• $E=s\cdot d\cdot \left(1-d^2+2c_1\right)+\left(1-d^2\right)\cdot \left(s+d\right)$

• $F=2s\cdot c_1\cdot d^2+\left(1-d^2\right)\cdot \left(s+2\right)$

Taking into account these assumptions the partial derivatives of Π₁ and V₂ take the following forms:

$\frac{\partial {\Pi }_1}{\partial p_1}=\frac{A-B\cdot p_1+C\cdot p_2}{{\left(1-d^2\right)}^2}$ and $\frac{\partial V_2}{\partial p_2}=\frac{D+E\cdot p_1-F\cdot p_2}{{\left(1-d^2\right)}^2}$

Both players follow the same strategy to decide their production quantities (homogeneous players) and they are characterized as bounded rational players. According to the existing literature it means that the first private company decides its price following a mechanism that is described by the equation:

$\frac{p_1\left(t+1\right)-p_1\left(t\right)}{p_1\left(t\right)}=k_1\cdot \frac{\partial {\Pi }_1}{\partial p_1}$ (14)

and the second player (semi-public company) who is also a bounded rational player following a similar mechanism that is given by the equation:

$\frac{p_2\left(t+1\right)-p_2\left(t\right)}{p_2\left(t\right)}=k_2\cdot \frac{\partial V_2}{\partial p_2}$ (15)

Through this mechanism each player increases his level of adaptation when his marginal utility is positive or decreases his level when his marginal utility is negative, where k is considered as the speed of adjustment for two players. The parameter k_i is positive (k_i > 0 with i = 1, 2) and gives the extend variation price of each i player is following a given signal. The dynamical system of the players is described by:

$\{\begin{array}{l}{p}_1\left(t+1\right)=p_1\left(t\right)+k_1\cdot p_1\left(t\right)\cdot \frac{\partial {\Pi }_1}{\partial p_1}\\ p_2\left(t+1\right)=p_2\left(t\right)+k_2\cdot p_2\left(t\right)\cdot \frac{\partial V_2}{\partial p_2}\end{array}$ (16)

The dynamical analysis of this system Equation (16) that is contained at the next sub-section is focusing on the parameters k₁ and k₂.

2.2. Dynamical Analysis

The dynamical analysis of the discrete dynamical system of Equation (16) involves the finding of equilibrium positions and the studying their stability. The ultimate goal of this algebraic study is to formulate a proposition that will be the stability condition of the Nash Equilibrium position. These algebraic results are verified and visualized doing some numerical simulations using the program of Mathematica.

2.2.1. The Equilibrium Positions

The equilibrium positions of the dynamical system Equation (16) are the nonnegative solutions of the following algebraic system:

$\{\begin{array}{l}{p}_1^{*}\cdot \frac{\partial {\Pi }_1}{\partial p_1}=0\\ p_2^{*}\cdot \frac{\partial V_2}{\partial p_2}=0\end{array}$ (17)

which is obtained by setting: $p_i\left(t+1\right)=p_i\left(t\right)=p_i^{*}$ for i = 1, 2.

• If $p_1^{*}=p_2^{*}=0$, then the equilibrium position is the point:

$E_0=\left(0,0\right)$ (18)

• If $p_1^{*}=0$ and $\frac{\partial V_2}{\partial p_2}=0$, then $p_2^{*}=\frac{D}{F}$ and the equilibrium position is the point:

$E_1=\left(0,\frac{D}{F}\right)$ (19)

• If $p_2^{*}=0$ and $\frac{\partial {\Pi }_1}{\partial p_1}=0$, then $p_1^{*}=\frac{A}{B}$ and the equilibrium position is the point:

$E_2=\left(\frac{A}{B},0\right)$ (20)

• If then the following system is obtained:

$\{\begin{array}{l}A-B\cdot p_1^{*}+C\cdot p_2^{*}=0\\ D+E\cdot p_1^{*}-F\cdot p_2^{*}=0\end{array}$ (21)

and the nonnegative solution of this algebraic gives the Nash Equilibrium position $E_{*}=\left(p_1^{*},p_2^{*}\right)$ where:

$p_1^{*}=\frac{A\cdot E+B\cdot D}{B\cdot F-C\cdot E}$ (22)

and:

$p_2^{*}=\frac{A}{B}+\frac{C\cdot \left(A\cdot E+B\cdot D\right)}{B\cdot \left(B\cdot F-C\cdot E\right)}$ (23)

From the three equilibrium positions E₁, E₂ and $E_{*}$, whereas A, B, C, D, E, F > 0 and $p_1^{*},p_2^{*}>0$, the following conditions arise:

$A\cdot E+B\cdot D>0$ (24)

$B\cdot F-C\cdot E>0$ (25)

2.2.2. Stability of Equilibrium Positions

To study the stability of the equilibrium positions the Jacobian matrix of the dynamical system Equation (16) is needed, which is the matrix:

$J\left(p_1^{*},p_2^{*}\right)=\left[\begin{array}{cc}{f}_{p_1}& f_{p_2}\\ g_{p_1}& g_{p_2}\end{array}\right]$ (26)

where:

$\begin{array}{l}f\left(p_1,p_2\right)=p_1+k_1\cdot p_1\cdot \frac{\partial {\Pi }_1}{\partial p_1}\\ g\left(p_1,p_2\right)=p_2+k_2\cdot p_2\cdot \frac{\partial V_2}{\partial p_2}\end{array}$ (27)

and as a result the Jacobian matrix of game’s discrete dynamical system Equation (16) is the following matrix:

$J\left(p_1^{*},p_2^{*}\right)=\left[\begin{array}{cc}1+k_1\cdot \left(\frac{\partial {\Pi }_1}{\partial p_1}+p_1^{*}\cdot \frac{{\partial }^2{\Pi }_1}{\partial p_1^2}\right)& k_1\cdot p_1^{*}\cdot \frac{{\partial }^2{\Pi }_1}{\partial p_1\partial p_2}\\ k_2\cdot p_2^{*}\cdot \frac{{\partial }^2{V}_2}{\partial p_2\partial p_1}& 1+k_2\cdot \left(\frac{\partial {\Pi }_2}{\partial p_2}+p_2^{*}\cdot \frac{{\partial }^2{V}_2}{\partial p_2^2}\right)\end{array}\right]$ (28)

For the E₀ the Jacobian matrix becomes as:

$J\left(E_0\right)=\left[\begin{array}{cc}1+k_1\cdot \frac{\partial {\Pi }_1}{\partial p_1}\left(E_0\right)& 0\\ 0& 1+k_2\cdot \frac{\partial V_2}{\partial p_2}\left(E_0\right)\end{array}\right]\underset{H=1+k_2\cdot \frac{\partial V_2}{\partial p_2}\left(E_0\right)}{\overset{G=1+k_1\cdot \frac{\partial {\Pi }_1}{\partial p_1}\left(E_0\right)}{=}}\left[\begin{array}{cc}G& 0\\ 0& H\end{array}\right]$ (29)

with Tr = G + H and Det = G × H.

The characteristic equation of J(E₀) is the following:

$r^2-TrJ\left(E_0\right)\cdot r+DetJ\left(E_0\right)=0$

$\iff r^2-\left(G+H\right)\cdot r+G\cdot H=0\iff \left(r-G\right)\cdot \left(r-H\right)=0$

The eigenvalues of J(E₀) are:

$r_1=G=1+k_1\cdot \frac{A}{{\left(1-d^2\right)}^2}$ and $r_2=H=1+k_2\cdot \frac{D}{{\left(1-d^2\right)}^2}$

As it is clearly seems that $\left|r_1\right|,\left|r_2\right|>1$, because A, D > 0. It means that the equilibrium position E₀ can be characterized as unstable.

For the position E₁ the Jacobian matrix becomes as:

$J\left(E_1\right)=\left[\begin{array}{cc}1+k_1\cdot \frac{\partial {\Pi }_1}{\partial p_1}\left(E_1\right)& 0\\ k_2\cdot p_2^{*}\cdot \frac{{\partial }^2{V}_2}{\partial p_2\partial p_1}\left(E_1\right)& 1+k_2\cdot p_2^{*}\cdot \frac{{\partial }^2{V}_2}{\partial p_2^2}\left(E_1\right)\end{array}\right]\underset{K=1+k_2\cdot p_2^{*}\cdot \frac{{\partial }^2{V}_2}{\partial p_2^2}}{\overset{I=1+k_1\cdot \frac{\partial {\Pi }_1}{\partial p_1}\left(E_1\right)}{=}}\left[\begin{array}{cc}I& 0\\ J& K\end{array}\right]$ (30)

with Tr = I + K and Det = I × K.

From the characteristic equation of J(E₁), the nonnegative eigenvalues are found as:

$r_1=I=1+k_1\cdot \frac{A+C\cdot p_2^{*}}{{\left(1-d^2\right)}^2}$ and $r_2=K=1-k_2\cdot \frac{F\cdot p_2^{*}}{1-d^2}$

Since, A, C > 0, it’s clearly seems that $\left|r_1\right|>1$ and the E₁ equilibrium position is unstable.

For the position E₂ the Jacobian matrix becomes as:

$J\left(E_2\right)=\left[\begin{array}{cc}1+k_1\cdot p_1^{*}\cdot \frac{{\partial }^2{\Pi }_1}{\partial p_1^2}\left(E_2\right)& k_1\cdot p_1^{*}\cdot \frac{{\partial }^2{\Pi }_1}{\partial p_1\partial p_2}\left(E_2\right)\\ 0& 1+k_2\cdot \frac{\partial V_2}{\partial p_2}\left(E_2\right)\end{array}\right]\underset{N=1+k_2\cdot \frac{\partial V_2}{\partial p_2}\left(E_2\right)}{\overset{L=1+k_1\cdot p_1^{*}\cdot \frac{{\partial }^2{\Pi }_1}{\partial p_1^2}\left(E_2\right)}{=}}\left[\begin{array}{cc}L& M\\ 0& N\end{array}\right]$ (31)

with Tr = L + N and Det = L × N.

From the characteristic equation of J(E₂), the nonnegative eigenvalues are found as:

$r_1=L=1-k_1\cdot \frac{B\cdot p_1^{*}}{{\left(1-d^2\right)}^2}$ and $r_2=N=1+k_2\cdot \frac{D+E\cdot p_1^{*}}{{\left(1-d^2\right)}^2}$

Since, D, E > 0, it’s clearly seems that $\left|r_2\right|>1$ and the E₂ equilibrium position is also unstable.

For the Nash equilibrium $E_{*}$ the Jacobian matrix becomes as:

$J\left(E_{*}\right)=\left[\begin{array}{cc}1+k_1\cdot p_1^{*}\cdot \frac{{\partial }^2{\Pi }_1}{\partial p_1^2}& k_1\cdot p_1^{*}\cdot \frac{{\partial }^2{\Pi }_1}{\partial p_1\partial p_2}\\ k_2\cdot p_2^{*}\cdot \frac{{\partial }^2{V}_2}{\partial p_2\partial p_1}& 1+k_2\cdot p_2^{*}\cdot \frac{{\partial }^2{V}_2}{\partial p_2^2}\end{array}\right]$ (32)

with:

$Tr=2+k_1\cdot p_1^{*}\cdot \frac{{\partial }^2{\Pi }_1}{\partial p_1^2}+k_2\cdot p_2^{*}\cdot \frac{{\partial }^2{V}_2}{\partial p_2^2}$ (33)

and:

$Det=\left(1+k_1\cdot p_1^{*}\cdot \frac{{\partial }^2{\Pi }_1}{\partial p_1^2}\right)\left(1+k_2\cdot p_2^{*}\cdot \frac{{\partial }^2{V}_2}{\partial p_2^2}\right)-k_1\cdot k_2\cdot p_1^{*}\cdot p_2^{*}\cdot \frac{{\partial }^2{\Pi }_1}{\partial p_1\partial p_2}\cdot \frac{{\partial }^2{V}_2}{\partial p_2\partial p_1}$ (34)

To study the stability of Nash equilibrium three conditions are used that the equilibrium position is locally asymptotically stable when they are satisfied simultaneously [34] [35]:

$\begin{array}{l}(\text{i})\hfill \\ (\text{ii})\hfill \\ (\text{iii})\hfill \end{array}\begin{array}{l}1-Det>0\hfill \\ 1-Tr+Det>0\hfill \\ 1+Tr+Det>0\hfill \end{array}$ (35)

It’s easy to find that the first inequality (i) becomes as:

$\begin{array}{l}1-Det>0\\ \iff k_1\cdot k_2\cdot p_1^{*}\cdot p_2^{*}\cdot \left(B\cdot F-C\cdot E\right)-\left(k_1\cdot p_1^{*}\cdot B+k_2\cdot p_2^{*}\cdot F\right)\cdot {\left(1-d^2\right)}^2<0\end{array}$ (36)

Also, the inequality (ii) gives:

$\begin{array}{l}1-Tr+Det>0\\ \iff k_1\cdot k_2\cdot p_1^{*}\cdot p_1^{*}\cdot \frac{B\cdot F-C\cdot E}{{\left(1-d^2\right)}^4}>0\end{array}$ (37)

and it’s always satisfied because of Equation (22).

Finally, the inequality (iii) becomes as:

$\begin{array}{l}1+Tr+Det>0\\ \iff k_1\cdot k_2\cdot p_1^{*}\cdot p_2^{*}\cdot \left(B\cdot F-C\cdot E\right)-2\left(k_1\cdot p_1^{*}\cdot B+k_2\cdot p_2^{*}\cdot F\right)\cdot {\left(1-d^2\right)}^2\\ +4{\left(1-d^2\right)}^4>0\end{array}$ (38)

Proposition:

The Nash equilibrium of the discrete dynamical system Equation (16) is locally asymptotically stable if:

$k_1\cdot k_2\cdot p_1^{*}\cdot p_2^{*}\cdot \left(B\cdot F-C\cdot E\right)-\left(k_1\cdot p_1^{*}\cdot B+k_2\cdot p_2^{*}\cdot F\right)\cdot {\left(1-d^2\right)}^2<0$

and:

$\begin{array}{l}{k}_1\cdot k_2\cdot p_1^{*}\cdot p_2^{*}\cdot \left(B\cdot F-C\cdot E\right)-2\left(k_1\cdot p_1^{*}\cdot B+k_2\cdot p_2^{*}\cdot F\right)\cdot {\left(1-d^2\right)}^2\\ +4{\left(1-d^2\right)}^4>0\end{array}$

where:

$p_1^{*}=\frac{A\cdot E+B\cdot D}{B\cdot F-C\cdot E}$ and $p_2^{*}=\frac{A}{B}+\frac{C\cdot \left(A\cdot E+B\cdot D\right)}{B\cdot \left(B\cdot F-C\cdot E\right)}$

3. Numerical Simulations

In this section the previous algebraic results are verified and the chaotic trajectories that appear when the parameters k₁ and k₂ take values outside the stability space are visualized. At first the stability region between the parameters k₁ and k₂ when all other game’s parameters take fixed values is plotted. Numerical simulations are made for the common speed of adjustment of two players (k₁ = k₂). Plotting the bifurcation diagrams of p₁ and p₂ with respect to the common parameter k (speed of adjustment), the stability interval is verified and the chaotic evolution through period doubling is highlighted. Strange attractors and Lyapunov numbers’ graphs reveal the chaotic behavior of the system when the parameter k takes larger values making the system unpredictable. Also, the sensitivity analysis in initial conditions [36] is made for these large values of the parameter k showing that only a small change on the coordinates of the system’s initial conditions can cause large differences in system’s behavior giving an evidence for chaos.

3.1. Stability Space between k₁ and k₂

The stability space between the two parameters k₁ and k₂ setting fixed values to the other parameters can created using the stability conditions of Proposition. These two inequalities Equation (36) and Equation (38) depend on the parameters α, c₁, c₂, c₄, d and s. Specific values are needed to be chosen for these parameters. For example the values of the parameters α = 5, c₁ = 1, c₂ = 0.5, c₄ = 1 and d = 0.5 are chosen. Taking into account that the percentage of the Greek State that participates in the Hellenic Petroleum Group is 35.5% it is necessary to set the value of 0.355 to the parameter s, which expresses this percentage. As a result the stability conditions take the following more specific forms:

$32.69\cdot k_1\cdot k_2-4.21\cdot k_1-3.18\cdot k_2<0$ (39)

$32.69\cdot k_1\cdot k_2-8.42\cdot k_1-6.36\cdot k_2+1.266>0$ (40)

Figure 3 shows the stability region between the parameters k₁ and k₂ as a result of Equation (39) and Equation (40).

Figure 3. Stability region between the parametersk₁ (horizontal axis) and k₂ (vertical axis) for α = 5, c₁ = 1, c₂ = 0.5, c₄ = 1, d = 0.5 and s = 0.355.

3.2. Focusing on the Common Speed of Adjustment k (k₁ = k₂)

At first, setting the previous specific values of the parameters α = 5, c₁ = 1, c₂ = 0.5,c₄ = 1, d = 0.5 and s = 0.355 the Nash Equilibrium becomes as:

$p_1^{*}\simeq 2.14$ and $p_2^{*}\simeq 2.91$ $\Rightarrow E_{*}\left(p_1^{*},p_2^{*}\right)\equiv E_{*}\left(2.14,2.91\right)$ (41)

and for the stability conditions it means that the parameter k must take values into the interval:

$k\in \left(0,0.097\right)$ (42)

This algebraic result is verified by the bifurcation diagrams of p₁ and p₂ with respect to the parameter k (Figure 4 and Figure 5). As it seems there is a locally asymptotically stable orbit until the value of 0.097 for the parameter k and for larger values of the parameter k the system’s behavior becomes chaotic and unpredictable.

(a)(b)

Figure 4. Bifurcation diagrams with respect to the parameter k against the variables p₁ (a) and p₂ (b) with 400 iterations of the map Equation (16) for α = 5, c₁ = 1, c₂ = 0.5, c₄ = 1, d = 0.5 ands = 0.355.

This chaotic trajectory can create strange attractors (Figure 6) for higher values of the parameter k like k = 0.14 (outside the stability space). Also, computing the Lyapunov numbers (Figure 7) for this value of the parameter k (k = 0.14) and setting the same fixed values for the other parameters, it seems that they are getting over the value of 1 as an evidence for the chaotic trajectory.

As it is known a characteristic property of chaotic trajectory is that the system becomes sensitive in initial conditions, which allows only a small change on the coordinates of initial conditions to change completely the system’s behavior. For example, choosing two different initial conditions (0.1, 0.1) and (0.101, 0.1) after a small change at the p₁—coordinate and plotting the time series of the system (Figure 8 and Figure 9) it seems that at the beginning the time series are indistinguishable, but after a number of iterations, the difference between them builds up rapidly.

Figure 5. The two bifurcation diagrams of Figure 4 in one.

Figure 6. Phase portrait (strange attractor) of the orbit of (0.1, 0.1) with 8000 iterations of the map Equation (16) for α = 5, c₁ = 1, c₂ = 0.5, c₄ = 1, d = 0.5, s = 0.355 and k = 0.14.

Figure 7. Lyapunov numbers of the orbit of (0.1, 0.1) with 8000 iterations of the map Equation (16) for α = 5, c₁ = 1, c₂ = 0.5, c₄ = 1, d = 0.5, s = 0.355 and k = 0.14.

(a)(b)

Figure 8. Sensitive dependence on initial conditions for p₁—coordinate plotted against the time: the orbit of (0.1, 0.1) (a) and the orbit of (0.101, 0.1) (b) of the system Equation (16) for α = 5, c₁ = 1, c₂ = 0.5, c₄ = 1, d = 0.5, s = 0.355 and k = 0.14.

Figure 9. The two bifurcation diagrams of Figure 8 in one.

4. Chaos Control Applying the d-Backtest Method

The d-Backtest is a parameters’ optimization method developed initially for financial applications [37] [38] [39], further more it has also been utilized in the research of the chaotic behavior of a duopoly, too [40]. In both cases the method attempts to identify backtesting periods that are optimal for parameter selection for a future period. for this application, the method utilizes various length verification periods and fixed length validation periods for the evaluation of the backtesting periods, using various methods that look into metrics such as profit factors, profits, extended profit margin and more [39]. This research uses the d-Backtest method for selecting values fork, dynamically for each player [40].

4.1. The Application of the d-Backtest Method

In this research the d-Backtest method is used for selecting values of the parameter k, dynamically and separately for each player [40]. Each player’s goal is to decrease the deviation in prices by choosing appropriate values for k_i in each step with the help of the d-Backtest method. In this case the system’s difference equations were different as the second player takes into account a percentage of the social welfare that it is containing in his utility function. Apart from the above differentiation, the use of d-Backtest remains the same with the metric being the root mean square deviation of the prices over a period:

$RMSD=\sqrt{\frac{{\displaystyle {\sum }_{i=1}^n{\left(x_i-\bar{x}\right)}^2}}{n-1}}$ (43)

This metric of variability is used for evaluating groups of backtesting periods by two families of methods, a simple average and an exponential average method. Each method family contains variations were the occurrences of backtesting periods are taken into account or the verification periods are longer than 1 and 2 steps or a combination of these two variations.

(a)(b)(c)(d)

Figure 10. (a) Price values with initial k = 0.147 for each player, (b) k values chosen by d-Backtest with initialk = 0.147 for each player, (c) Backtesting periods selecting by the d-Backtest method with initial k = 0.147 for each player, (d) Phase diagram of price values with initial k = 0.147 for each player.

(a)(b)(c)(d)

Figure 11. (a) Price values with initial k = 0.3 for each player, (b) k values chosen by d-Backtest with initial k = 0.3 for each player, (c) Backtesting periods selecting by the d-Backtest method with initialk = 0.3 for each player, (d) Phase diagram of price values with initial k = 0.3 for each player.

(a)(b)(c)(d)

Figure 12. (a) Price values with initialk = 0.7 for each player, (b) k values chosen by d-Backtest with initial k = 0.7 for each player, (c) Backtesting periods selecting by the d-Backtest method with initial k = 0.7 for each player, (d) Phase diagram of price values with initial k = 0.7 for each player.

Each player calculates the backtesting prices separately from each other, using the others’ past prices in the formulas calculating their prices over past periods. The prices have to be bound to a lowest value of 0.001 so that no negative prices appear. For each step each player calculates the prices for all 30 backtesting periods and all possible k_i values in [0.001, 1.0] with a step of 0.001. Then, they compute the variability metric for each period and k_i value and provide the data to the d-Backtesting method that will choose the right method for sorting the baktesting periods, then the right baktesting period and then the right k_i value for the next step. Each player starts the game with an initial price of 0.1 and a value of 0.1 for the parameter k for the first 10 steps, until there are enough data for the d-Backtest method to work with after that.

4.2. The d-Backtest Method’s Results

Games for initial values of k as 0.147 (Figure 10), 0.3 (Figure 11) and 0.7 (Figure 12) and 500 steps for each one were ran. The results are shown in price values, k values and backtesting periods through the following figures:

The results resemble to those of Sarafopoulos et al. [40]. In every game the price manage to get closer to the theoretical equilibrium regardless of the initial values of k. this also happens regardless of the fact that the values chosen for k are not always in the stability zone or most of the time not in stability zone. The k values for the second player are more diverse than those of the first player, with the value of 0.5 appearing often as the d-Backtest method chooses the middle value when there is not different in the value space examined. The backtesting periods are mostly chosen with values of 1 or 2 which happens because the calculated variability is usually lower in small backtesting periods.

5. Conclusion

In this paper, the dynamics of a differentiated Bertrand duopoly with homogeneous expectations, linear demand and asymmetric cost functions (regression analysis) are analyzed. The different cases, in which one of the players is not a completely private company and care about a percentage of the social welfare due to the participation of the public capital into it, are studied. By assuming that at each time period, each firm maximizes its expected utility function Π₁ (for the first player) or V₂ (for the second player) under bounded rationality expectation, a discrete dynamic system was obtained. The existence and stability of equilibrium of this system are studied. It is numerically shown that the model gives chaotic and unpredictable trajectories. The main result is that higher values of the speed of adjustment may destabilize the Bertrand-Nash equilibrium. Finally, in cases where the players choose these values of the parameter k for which, as algebraically and graphically proved, the discrete dynamical system of the Bertrand-type model behaves chaotically, the d-Backtest method was applied, giving dynamic values for the k parameter to each player for each time period, through which the system returns in locally asymptotically stable Nash Equilibrium.

Acknowledgements

The authors thank the anonymous reviewers for interesting remarks and useful suggestions.

Funding

This research is carried out/funded in the context of the project “Nonlinear Dynamics in an Oligopoly Market: Chaotic Behavior, Complexity Control Using Machine Learning and Application in an Oligopoly of the Greek Market” (MIS 5049905) under the call for proposals “researchers’ support with an emphasis on young researchers—2^nd cycle”. The project is co-financed by Greece and the European Union (European Social Fund—ESF) by the Operational Programme “Human Resources Development, Education and Lifelong Learning 2014-2020”.

NOTES

¹The coefficients of determination of the regressions are relatively small 0.51 and 0.39 respectively. This is due to the fact that the cost functions of the companies depend on many unpublished factors.

Conflicts of Interest

The authors declare no conflicts of interest regarding the publication of this paper.

References

[1]	Agiza, H.N. (1999) On the Analysis of Stability, Bifurcation, Chaos and Chaos Control of Kopel Map. Chaos, Solitons & Fractals, 10, 1909-1916. https://doi.org/10.1016/S0960-0779(98)00210-0
[2]	Agiza, H.N., Hegazi, A.S. and Elsadany, A.A. (2002) Complex Dynamics and Synchronization of Duopoly Game with Bounded Rationality. Mathematics and Computers in Simulation, 58, 133-146. https://doi.org/10.1016/S0378-4754(01)00347-0
[3]	Agliar, A., Gardini, L. and Puu, T. (2005) Some Global Bifurcations Related to the Appearance of Closed Invariant Curves. Mathematics and Computers in Simulation, 68, 201-219. https://doi.org/10.1016/j.matcom.2004.12.003
[4]	Agliari, A., Gardini, L. and Puu, T. (2006) Global Bifurcations in Duopoly When the Cournot Point Is Destabilized via a Subcritical Neimark Bifurcation. International Game Theory Review, 8, 1-20. https://doi.org/10.1142/S0219198906000758
[5]	Bishi, G.I. and Kopel, M. (2001) Equilibrium Selection in a Nonlinear Duopoly Game with Adaptive Expectations. Journal of Economic Behavior & Organization, 46, 73-100. https://doi.org/10.1016/S0167-2681(01)00188-3
[6]	Kopel, M. (1996) Simple and Complex Adjustment Dynamics in Cournot Duopoly Models. Chaos, Solitons & Fractals, 7, 2031-2048. https://doi.org/10.1016/S0960-0779(96)00070-7
[7]	Puu, T. (1998) The Chaotic Duopolists Revisited. Journal of Economic Beahavior & Organization, 33, 385-394. https://doi.org/10.1016/S0167-2681(97)00064-4
[8]	Puu, T. (2005) Complex Oligopoly Dynamics. In: Lines, M., Ed., Nonlinear Dynamical Systems in Economics, Springer, New York, 165-186. https://doi.org/10.1007/3-211-38043-4_6
[9]	Sarafopoulos, G. (2015) On the Dynamics of a Duopoly Game with Differentiated Goods. Procedia Economics and Finance, 19, 146-153. https://doi.org/10.1016/S2212-5671(15)00016-7
[10]	Sarafopoulos, G. (2015) Complexity in a Duopoly Game with Homogeneous Players, Convex, Log-Linear Demand and Quadratic Cost Functions. Procedia Economics and Finance, 33, 358-366. https://doi.org/10.1016/S2212-5671(15)01720-7
[11]	Zhang, J., Da, Q. and Wang, Y. (2009) The Dynamics of Bertrand Model with Bounded Rationality. Chaos, Solitons and Fractals, 39, 2048-2055. https://doi.org/10.1016/j.chaos.2007.06.056
[12]	Agiza, H.N. and Elsadany, A.A. (2003) Nonlinear Dynamics in the Cournot Duopoly Game with Heterogeneous Players. Physica A, 320, 512-524. https://doi.org/10.1016/S0378-4371(02)01648-5
[13]	Agiza, H.N. and Elsadany, A.A. (2004) Chaotic Dynamics in Nonlinear Duopoly Game with Heterogenous Players. Applied Mathematics and Computation, 149, 843-860. https://doi.org/10.1016/S0096-3003(03)00190-5
[14]	Haan, W.J.D. (2001) The Importance of the Number of Different Agents in a Heterogeneous Asset—Pricing Model. Journal of Economic Dynamic and Control, 25, 721-746. https://doi.org/10.1016/S0165-1889(00)00038-5
[15]	Hommes, C.H. (2006) Heterogeneous Agent Models in Economics and Finance. In: Tesfatsion, L. and Judd, K.L., Eds., Handbook of Computational Economics, Agent-Based Computational Economics, Vol. 2, Elsevier Science B.V., Amsterdam, 1109-1186. https://doi.org/10.1016/S1574-0021(05)02023-X
[16]	Fanti, L. and Gori, L. (2012) The Dynamics of a Differentiated Duopoly with Quantity Competition. Economic Modelling, 29, 421-427. https://doi.org/10.1016/j.econmod.2011.11.010
[17]	Gao, Y. (2009) Complex Dynamics in Two Dimensional Noninvertible Map. Chaos Solitons & Fractals, 39, 1798-1810. https://doi.org/10.1016/j.chaos.2007.06.051
[18]	Sarafopoulos, G. and Papadopoulos, K. (2017) On a Cournot Duopoly Game with Differentiated Goods, Heterogeneous Expectations and a Cost Function Including Emission Costs. Scientific Bulletin—Economic Science, 16, 11-22. https://econpapers.repec.org/article/ptsjournl/y_3a2017_3ai_3a1_3ap_3a11-22.htm
[19]	Sarafopoulos, G. and Papadopoulos, K. (2019) Complexity in a Bertrand Duopoly Game with Heterogeneous Players and Differentiated Goods. In: Sykianakis, N., Polychronidou, P. and Karasavvoglou, A., Eds., Economic and Financial Challenges for Eastern Europe, Springer, Berlin, 15-26. https://doi.org/10.1007/978-3-030-12169-3_2
[20]	Sarafopoulos, G. and Papadopoulos, K. (2020) On a Bertrand Dynamic Game with Differentiated Goods, Heterogeneous Expectations and Asymmetric Cost Functions. In: Janowicz-Lomott, M., Łyskawa, K. and Polychronidou, P, Eds., Economic and Financial Challenges for Balkan and Eastern European Countries, Springer, Berlin, 223-241. https://doi.org/10.1007/978-3-030-39927-6_14
[21]	Tramontana, F. (2010) Heterogeneous Duopoly with Isoelastic Demand Function. Economic Modelling, 27, 350-357. https://doi.org/10.1016/j.econmod.2009.09.014
[22]	Zhang, J., Da, Q. and Wang, Y. (2007) Analysis of Nonlinear Duopoly Game with Heterogeneous Players. Economic Modelling, 24, 138-148. https://doi.org/10.1016/j.econmod.2006.06.007
[23]	Wu, W., Chen, Z. and Ip, W.H. (2010) Complex Nonlinear Dynamics and Controlling Chaos in a Cournot Duopoly Economic Model. Nonlinear Analysis: Real World Applications, 11, 4363-4377. https://doi.org/10.1016/j.nonrwa.2010.05.022
[24]	Baumol, W.J. and Quandt, R.E. (1964) Rules of Thumb and Optimally Imperfect Decisions. American Economic Review, 54, 23-46. https://www.jstor.org/stable/1810896
[25]	Singh, N. and Vives, X. (1984) Price and Quantity Competition in a Differentiated Duopoly. The RAND Journal of Economics, 15, 546-554. https://doi.org/10.2307/2555525
[26]	Puu, T. (1991) Chaos in Duopoly Pricing. Chaos, Solitons & Fractals, 1, 573-581. https://doi.org/10.1016/0960-0779(91)90045-B
[27]	Puu, T. (1995) The Chaotic Monopolist. Chaos, Solitons & Fractals, 5, 35-44. https://doi.org/10.1016/0960-0779(94)00206-6
[28]	Westerhoff, F. (2006) Nonlinear Expectation Formation, Endogenous Business Cycles and Stylized Facts. Studies in Nonlinear Dynamics and Econometrics, 10, Article No. 4. https://doi.org/10.2202/1558-3708.1324
[29]	Naimzada, A.K. and Ricchiuti, G. (2008) Complex Dynamics in a Monopoly with a Rule of Thumb. Applied Mathematics and Computation, 203, 921-925. https://doi.org/10.1016/j.amc.2008.04.020
[30]	Asksar, S.S. (2013) On Complex Dynamics of Monopoly Market. Economic Modelling, 31, 586-589. https://doi.org/10.1016/j.econmod.2012.12.025
[31]	Askar, S.S. (2014) Complex Dynamic Properties of Cournot Duopoly Games with Convex and Log-Concave Demand Function. Operations Research Letters, 42, 85-90. https://doi.org/10.1016/j.orl.2013.12.006
[32]	Elsadany, A. and Awad, A.M. (2016) Nonlinear Dynamics of Cournot Duopoly Game with Social Welfare. Electronic Journal of Mathematician Analysis and Applications, 4, 173-191. http://math-frac.org/Journals/EJMAA/Vol4(2)_July_2016/Vol4(2)_Papers/16_EJMAA_Vol4(2)_July_2016_pp_173-191.pdf
[33]	Naimzada, A.K. and Sbragia, L. (2006) Oilgopoly Games with Nonlinear Demand and Cost Functions: Two Boundedly Rational Adjustment Processes. Chaos, Solitons & Fractals, 29, 707-722. https://doi.org/10.1016/j.chaos.2005.08.103
[34]	Gandolfo, G. (1997) Economic Dynamics. Springer, Berlin. https://www.springer.com/gp/book/9783642038624
[35]	Elaydi, S. (2005) An Introduction to Difference Equations. 3rd Edition, Springer-Verlag, New York.
[36]	Kulenovic, M. and Merino, O. (2002) Discrete Dynamical Systems and Difference Equations with Mathematica. Chapman & Hall/CRC, London. https://doi.org/10.1201/9781420035353
[37]	Vezeris, D., Schinas, C. and Papaschinopoulos, G. (2018) Profitability Edge by Dynamic Back Testing Optimal Period Selection Technical Parameters Optimization in Trading Systems with Forecasting. Computational Economics, 51, 761-807. https://doi.org/10.1007/s10614-016-9640-x
[38]	Vezeris, D., Schinas, C., Kyrgos, T., Bizergianidou, V. and Karkanis, I. (2019) Optimizations of Backtesting Techniques in Automated High Frequency Trading Systems Using the d-Backtest PS Method. Computational Economics, 56, 975-1054. https://doi.org/10.1007/s10614-019-09956-1
[39]	Vezeris, D., Kyrgos, T., Karkanis, I. and Bizergianidou, V. (2020) Automated Trading Systems’ Evaluation Using d-Backtest PS Method and WM Ranking in Financial Markets. Investment Management and Financial Innovations, 17, 198-215. https://doi.org/10.21511/imfi.17(2).2020.16
[40]	Sarafopoulos, G., Drimpetas, E., Papadopoulos, K. and Vezeris, D. (2021) Chaotic Behavior in Duopoly Market and Application of the d-Backtest Method. 14th Chaotic Modeling and Simulation International Conference CHAOS 2021, Athens, 8-11 June 2021.