A Continuous Model of Marital Relations with Stochastic Differential Equations

Abstract

Marital relations depend on many factors which can increase the amount of satisfaction or unhappiness in the relation. A large percentage of marriages end up in divorce. While there are many studies about the causes of divorce and how to prevent it, there are very few mathematical models dealing with marital relations. In this paper, we present a continuous model based on the ideas presented by Gottman and coauthors. We show that the type of influence functions that describe the interaction between husband and wife is critical in determining the outcome of a marriage. We also introduce stochasticity into the model to account for the many factors that affect the marriage and that are not easily quantified, such as economic climate, work stress, and family relations. We show that these factors are able to change the equilibrium state of the couple.

1. Introduction

Marital relations are very complicated. The satisfaction or lack of satisfaction in a couple depends on the personality of each partner and on the interaction between them, but also on many other external and internal factors, such as the economic and health situations. Mathematical models together with empirical studies can help better understand the relations and determine which factors can improve the satisfaction the most. The first mathematical models dealing with the relation of a couple were linear models, such as the Romeo and Juliet model [1], given by two linear ordinary differential equations in which the rate of change of the satisfaction of each Romeo and Juliet is a linear combination of the satisfaction of both. This model was also studied in Reference [2] and references therein. Gosh [3] improved on this model by adding a as a third variable, the effort put by the male to change the relation. In Reference [4], an optimal control problem for the satisfaction is studied. In Reference [5], the appearance of chaos is studied in a similar model, and, in Reference [6], periodic solutions are analyzed.

In the work by Gottman et al. [7], they introduced a set of discrete equations in discrete time t to model the scores of a measured interaction between a married couple. A negativity/positivity score, called Gottman-Levenson score, is given to each turn at speech. The researchers categorized the affects displayed during the interaction using five positive codes (interest, validation, affection, humor, joy), ten negative affect codes (disgust, contempt, belligerence, domineering, anger, fear/tension, defensiveness, whining, sadness, stonewalling), and a neutral affect code. Positive affects had an associated positive integer, negative affects an associated negative integer, and the neutral affect had value zero. The sum of the affects, $W_t$ for the wife and $H_t$ for the husband, give a rough measure of the satisfaction during the tth turn at speech. The more positive the value, the more satisfied the person will be. It was assumed that each person’s score is determined solely by his/her own and his/her partner’s immediately previous score and thus the equations have the form:

$\begin{array}{ccc}\hfill W_{t+1}& =& f(W_t,H_t)\hfill \\ \hfill H_{t+1}& =& g(W_{t+1},H_t).\hfill \end{array}$

The functions f and g have to be determined and, in general, are different for each couple. Some simplifying assumptions are made on these functions. They assume that the past two scores contribute separately and that the effects can be added together. Hence, a person’s score is regarded as the sum of two components: (1) interpersonal influence from spouse to spouse, termed influenced behavior, and (2) the individual’s own dynamics, termed uninfluenced behavior. Uninfluenced behavior depends on each individual’s previous score only, and influenced behavior depends on the score for the partner’s last turn of speech.

The uninfluenced part has the form:

$P_{t+1}={\tilde{r}}_P{P}_t+a_P,$

where $P=H$ for the husband, and $P=W$ for the wife. P is Gottman-Levenson variable score at term at speech t; $r_P$, the emotional inertia, is the rate at which the individual returns to the uninfluenced equilibrium point, and $a_P$ is a constant specific to the individual. If we call P the uninfluenced equilibrium (fixed) point, then P and $a_p$ are related by

$P=a_p/(1-r_p).$

They introduce the effect of the husband on the wife, and vice-versa, by defining influence functions, $I_{WH}$ for the influence of the wife on the husband and $I_{HW}$ for the one of the husband on the wife. So, for example, their equation for the satisfaction of the husband is

(1)$H_{t+1}={\tilde{r}}_H{H}_t+a_H+I_{WH}\left(W_t\right).$

They consider several different types of personalities for both the husband and the wife: volatile, validating, conflict-avoiding, hostile, and hostile-detached. In all cases, they consider two types of influence functions, bilinear and ojive or piece-wise constant. The bilinear influence functions have the form:

(2)$I\left(P\right)=\left\{\begin{array}{cc}{m}_{m}P,\hfill & P<0\hfill \\ m_{p}P,\hfill & 0\le P\hfill \end{array}\right.,$

where $m_m$ and $m_p$ are the corresponding slopes. The ojive piece-wise constant influence functions have the form:

(3)$I\left(P\right)=\left\{\begin{array}{cc}{I}_{-2},\hfill & P<T_{-}\hfill \\ I_{-1},\hfill & T_{-}\le P<0\hfill \\ I_{+1},\hfill & 0\le P<T_{+}\hfill \\ I_{+2},\hfill & T_{+}\le P\hfill \end{array}\right.,$

where $I_{-2},I_{-1},I_{+1}$, and $I_{+2}$ are constants.

We developed a modification of Gottman’s work which consists of a couple of differential equations with more general influence functions and a stochastic term which models possible external influences. We first analyze equilibrium points and their stability for the deterministic continuous case and then study a stochastic version which uses white noise. We study only the case where both husband and wife are validating, which means that P and $I\left(P\right)$ have the same sign. We vary the influence functions and external contributions.

1.1. Continuous Model

Deterministic Case

We consider a continuous time model based on the model proposed by Gottman et al. [7].

(4)$\begin{array}{cc}\hfill \frac{dH}{dt}=& a_H+r_{H}H+I_{WH}\left(W\right)\hfill \\ \hfill \frac{dW}{dt}=& a_W+r_{W}W+I_{HW}\left(H\right).\hfill \end{array}$

We fix the values of $r_H,r_W,a_H$, and $a_W$ so that the uninfluenced equilibrium states coincide with the ones in Table 10.1 of Reference [7]. The discrete model in Reference [7] is Euler’s approximation of (4) with time step equal to 1, as in Reference [8]. For our first equation,

(5)$\frac{dH}{dt}\simeq \frac{H_{t+1}-H_t}{1}=a_H+r_H{H}_t+I_{WH}\left(W_t\right).$

Hence,

(6)$H_{t+1}=a_H+(r_H+1)H_t+I_{WH}\left(W_t\right).$

Therefore, comparing with Equation (1), $r_H={\tilde{r}}_H-1$, and $r_W={\tilde{r}}_W-1$. The parameters $a_H$ and $a_W$ depend on each person, and they are the same for the discrete and the continuous models. They are connected to the uninfluenced steady states in the following way. If there is no influence,

(7)$\frac{dW}{dt}=a_W+r_{W}W;$

hence, the uninfluenced equilibrium point for the wife is given by

(8)$W=\frac{-a_W}{r_W}.$

We study the effect of four types of influence functions. We will consider only the case when both husband and wife are validating. According to Gottman: “Validators have an influence function that creates an influence towards negativity in a spouse if the partner’s behavior is negative and an influence towards positivity if the partner’s behavior is positive.” In Figure 1, we show an example with the values we will be using in the numerical simulations, which were taken from Reference [7].

We will use as example the values of the inertia terms for husband and wife given in Gottman’s Table 10.1, which are ${\tilde{r}}_H=0.37$ and ${\tilde{r}}_W=0.14$. Hence, $r_H=-0.63$ and $r_W=-0.86$, both negative. All the other cases that are considered in Reference [7] can be dealt with in a similar way.

The first influence function considered is Gottman’s bilinear function (2):

(9)$I_{WH}\left(W\right)=\left\{\begin{array}{cc}{M}_{-}W,\hfill & W<0\hfill \\ M_{+}W,\hfill & 0\le W\hfill \end{array}\right.,$

(10)$I_{HW}\left(H\right)=\left\{\begin{array}{cc}{m}_{-}H,\hfill & H<0\hfill \\ m_{+}H,\hfill & 0\le H\hfill \end{array}\right..$

Since it has been observed [7] that, for validating couples, when one spouse is unhappy, his or her effect on the other spouse is larger than when he or she is happy, and the values of the slopes for the negative parts are larger than for the positive parts.

The second one consists of piece-wise constant functions, just like the ojive function of Reference [7]:

(11)$I_{WH}\left(W\right)=\left\{\begin{array}{cc}{I}_{-2}^W,\hfill & W<W_{-2}\hfill \\ I_{-1}^W,\hfill & W_{-2}\le W<W_{-1}\hfill \\ 0,\hfill & W-1\le W<W_1\hfill \\ I_1^W,\hfill & W_1\le W<W_2\hfill \\ I_2^W,\hfill & W_2\le W\hfill \end{array}\right.,$

(12)$I_{HW}\left(H\right)=\left\{\begin{array}{cc}{I}_{-2}^H,\hfill & H<H_{-2}\hfill \\ I_{-1}^H,\hfill & H_{-2}\le H<H_{-1}\hfill \\ 0,\hfill & H_{-1}\le H<H_1\hfill \\ I_1^H,\hfill & H_1\le H<H_2\hfill \\ I_2^H,\hfill & H_2\le H\hfill \end{array}\right..$

$I_j^W$ are constants that define the influence function of the wife on the husband, and $I_j^H$ are constants that define the influence function of the husband on the wife. In the analysis done in Reference [7], these functions are non-decreasing, and $H_{-1}=H_1=0$.

In this work, we consider other scenarios based on conversations with a family therapist and on specific cases studied by Gottman, where $I\left(P\right)$ is better fitted with more general functions and not necessarily bilinear or piece-wise constant. For instance, we assume $H_{-1}<0<H_1$, which means that small values of the happiness of one spouse do not affect the other spouse’s happiness. We also study the possibility that, for larger values of the satisfaction (or dissatisfaction), there may be a saturation effect, and the influence function decreases or even goes to zero. The corresponding nullclines can be seen in Figure 2, Figure 3 and Figure 4.

We introduce a third type of influence function where we assume that, for small values of the partner’s satisfaction, the person’s satisfaction increases with a linear rate up to a point and then the rate decreases linearly until, for a certain value, the rate is zero, and similarly for negative values of the satisfaction. In this way, we introduce a saturation of the influence of one spouse on the other. The influence functions are:

(13)$I_{WH}\left(W\right)=\left\{\begin{array}{cc}0,\hfill & W<W_{-2}\hfill \\ \frac{-W_{-2}{I}_{-1}^W}{W_{-1}-W_{-2}}+\frac{I_{-1}^W}{W_{-1}-W_{-2}}W,\hfill & W_{-2}\le W<W_{-1}\hfill \\ \frac{I_{-1}^W}{W_{-1}}W,\hfill & W_{-1}\le W<0\hfill \\ \frac{I_1^W}{W_1}W,\hfill & 0\le W<W_1\hfill \\ \frac{-W_2{I}_{+1}^W}{W_1-W_2}+\frac{I_1^W}{W_1-W_2}W,\hfill & W_1\le W<W_2\hfill \\ 0,\hfill & W_2\le W\hfill \end{array}\right.$

(14)$I_{HW}\left(H\right)=\left\{\begin{array}{cc}0,\hfill & H<H_{-2}\hfill \\ \frac{-H_{-2}{I}_{-1}^H}{H_{-1}-H_{-2}}+\frac{I_{-1}^H}{H_{-1}-H_{-2}}H,\hfill & H_{-2}\le H<H_{-1}\hfill \\ \frac{I_{-1}^H}{H_{-1}}H,\hfill & H_{-1}\le H<0\hfill \\ \frac{I_1^H}{H_1}H,\hfill & 0\le H<H_1\hfill \\ \frac{-H_2{I}_1^H}{H_1-H_2}+\frac{I_1^H}{H_1-H_2}H,\hfill & H_1\le H<H_2\hfill \\ 0,\hfill & H_2\le W\hfill \end{array}\right..$

The final type of influence functions that we consider are similar to those of the third type only replacing the piece-wise linear functions by piece-wise cubic functions so that the influence functions will not only be continuous but also smooth.

(15)$I_{WH}\left(W\right)=\left\{\begin{array}{cc}0,\hfill & W<W_{-2}\hfill \\ \frac{2}{3(W_{-2}-W_{-1})}{W}^3+\frac{3I_{-1}^W+2W_{-2}^2-4W_{-1}{W}_{-2}+2W_{-1}^2}{3W_{-2}^2-6W_{-1}{W}_{-2}+3W_{-1}^2}{W}^2,\hfill & W_{-2}\le W<W_{-1}\hfill \\ -\frac{2}{3W_{-1}}{W}^3+\frac{3I_{-1}^W+2W_{-1}^2}{3W_{-1}^2}{W}^2,\hfill & W_{-1}\le W<0\hfill \\ -\frac{2}{3W_1}{W}^3+\frac{3I_1^W+2W_1^2}{3W_1^2}{W}^2,\hfill & 0\le W<W_1\hfill \\ \frac{2}{3(W_2-W_1)}{W}^3+\frac{3I_1^W+2W_2^2-4W_1{W}_2+2W_1^2}{3W_2^2-6W_1{W}_2+3W_1^2}{W}^2,\hfill & W_1\le W<W_2\hfill \\ 0,\hfill & W_{+2}\le W\hfill \end{array}\right.,$

(16)$I_{HW}\left(H\right)=\left\{\begin{array}{cc}0,\hfill & H<H_{-2}\hfill \\ \frac{2}{3(H_{-2}-H_{-1})}{H}^3+\frac{3I_{-1}^H+2H_{-2}^2-4H_{-1}{H}_{-2}+2H_{-1}^2}{3H_{-2}^2-6H_{-1}{H}_{-2}+3H_{-1}^2}{H}^2,\hfill & H_{-2}\le H<H_{-1}\hfill \\ -\frac{2}{3H_{-1}}{H}^3+\frac{3I_{-1}^H+2H_{-1}^2}{3H_{-1}^2}{H}^2,\hfill & H_{-1}\le H<0\hfill \\ -\frac{2}{3H_1}{H}^3+\frac{3I_1^H+2H_1^2}{3H_1^2}{H}^2,\hfill & 0\le H<H_1\hfill \\ \frac{2}{3(H_2-H_1)}{H}^3+\frac{3I_1^H+2H_2^2-4H_1{H}_2+2H_1^2}{3H_2^2-6H_1{H}_2+3H_1^2}{H}^2,\hfill & H_1\le H<H_2\hfill \\ 0,\hfill & H_2\le H\hfill \end{array}\right..$

Graphs for these types of influence functions for fixed values of the parameters are given in Section 3. The values of the parameters we chose for all the new influence functions are only academic; they are not based on actual experiments with couples.

1.2. Stochastic Model

The satisfaction or happiness of an individual depends not only on the personality of the person and his or her spouse but also on economic factors, weather, family, political situation, and many other factors. These factors are very hard to quantify and have a large degree of variability and uncertainty. One way to include them in mathematical models is by using Ito Stochastic Differential Equations [9,10].

A class of frequently used stochastic processes is the Brownian one.

We will add a white noise term to our deterministic models.

(17)$\begin{array}{cc}\hfill dH=& (a_H+r_{H}H+I_{WH}\left(W\right))dt+s_{1}d{W}_H\left(t\right)\hfill \\ \hfill dW=& (a_W+r_{W}W+I_{HW}\left(H\right))dt+s_{2}d{W}_W\left(t\right).\hfill \end{array}$

Here, $W_H\left(t\right)$ and $W_W\left(t\right)$ are Brownian motions also known in Mathematics as Wiener Processes. Since we assume that the noise does not depend on either spouse’s satisfaction, we are adding the white noise multiplied by a constant $s_1$ or $s_2$.

2. Analysis of the Deterministic Case

In this section, we study the steady or equilibrium solutions of Equation (4) and their local stability for the four types of influence functions considered. The equilibrium solutions are found by determining the intersection of the nullclines of (4). The local stability of an equilibrium point is found by linearizing (4) around the equilibrium point and finding the eigenvalues of the Jacobian matrix. If all the eigenvalues have negative real part, then the equilibrium point is locally stable. Equilibrium points represent the state of happiness (if the equilibrium point is positive) or unhappiness (if the equilibrium point is negative) the married couple tends to approach.

Qualitative Analysis for the Bilinear Influence Function

We consider system (4) with influence function of bilinear type given as:

(18)$I_{WH}\left(W\right)=\left\{\begin{array}{cc}{M}_{-}W,\hfill & W<0\hfill \\ M_{+}W,\hfill & 0\le W\hfill \end{array}\right.,$

(19)$I_{HW}\left(H\right)=\left\{\begin{array}{cc}{m}_{-}H,\hfill & H<0\hfill \\ m_{+}H,\hfill & 0\le H\hfill \end{array}\right..$

These functions define different systems, depending on the quadrant.

(20)$\frac{dH}{dt}=\left\{\begin{array}{cc}{a}_H+r_{H}H+M_{-}W,\hfill & W<0\hfill \\ a_H+r_{H}H+M_{+}W,\hfill & 0\le W\hfill \end{array}\right.,$

(21)$\frac{dW}{dt}=\left\{\begin{array}{cc}{a}_W+r_{W}W+m_{-}H,\hfill & H<0\hfill \\ a_W+r_{W}W+m_{+}H,\hfill & 0\le H\hfill \end{array}\right..$

The equilibrium points for each case, $(H^{*},W^{*})$, can be written as:

(22)$H^{*}=\frac{a_{W}M-a_H{r}_W}{r_H{r}_W-mM},W^{*}=\frac{a_{H}m-a_W{r}_H}{r_H{r}_W-mM},$

where the values of m and M are taken, according to the quadrant under consideration, to be $m_{-}$ or $m_{+}$ and $M_{-}$ or $M_{+}$.

If the value of $W\left(orH\right)$ is positive, it means the wife (or husband) is happy. If $W\left(orH\right)$ is negative, it means the wife (or husband) is unhappy. Therefore, if the equilibrium point is in the first quadrant, then both wife and husband are happy. In the third quadrant, both are unhappy. In the second quadrant, the husband is unhappy, and the wife is happy; and, in the fourth quadrant, the wife is unhappy, and the husband is happy.

As usual, the stability of the equilibria is analyzed by linearizing the system around the equilibrium points. The corresponding Jacobians of the linearized system are:

$\left(\begin{array}{cc}{r}_H\hfill & M\hfill \\ m\hfill & r_W\hfill \end{array}\right).$

Therefore, $trace=r_H+r_W=p$, which is negative for the given r’s, $\mathrm{determinant}=r_W{r}_H-mM=q$; then, the discriminant is $\Delta =p^2-4q={(r_H-r_W)}^2+4mM$, and the eigenvalues are

${\lambda }_{\pm }=\frac{p\pm \sqrt{\Delta }}{2}.$

Therefore, the critical point of the linearized system will be:

• (1). A saddle point if:

  $\Delta >0,q<0\iff {(r_H-r_W)}^2+4mM>0\mathrm{and}{r}_H{r}_W<mM.$

• (2). A node if:

  $\Delta \ge 0,q>0\iff {(r_H-r_W)}^2+4mM\ge 0\mathrm{and}{r}_H{r}_W>mM.$

  It will be stable if $p<0$ or unstable if $p>0$.

• (3). A spiral if:

  $\Delta <0,p\ne 0\iff {(r_H-r_W)}^2+4mM<0\mathrm{and}{r}_H+r_W\ne 0;$

  it will be stable if $p<0$ and unstable if $p>0$.

• (4). A center if:

  $\Delta <0,p=0\iff {(r_H-r_W)}^2+4mM<0\mathrm{and}{r}_H+r_W=0.$

• (5). A proper or improper node if:

  $\Delta =0\iff {(r_H-r_W)}^2+4mM=0.$

It is known (Hartman-Grobman Theorem) that the fixed points of the non-linear system and of its linearization will be of the same kind except when the eigenvalues have a real part equal to zero.

For the problem under investigation, we have $p<0$, and, for the validating model, $m>0,M>0$ always; therefore, $\Delta \ge 0$ and, using the above results, we have proved the following:

Theorem 1.

The equilibrium points, or steady states of system (4) for the bilinear influence function, are given by

$H^{*}=\frac{a_{W}M-a_H{r}_W}{r_H{r}_W-mM}{W}^{*}=\frac{a_{H}m-a_W{r}_H}{r_H{r}_W-mM},$

for the appropriate values of m and M on the corresponding quadrant. These points are:

• (i) 

  Saddle points (and therefore unstable) if and only if $r_H{r}_W<mM$, and ${(r_H-r_W)}^2+4mM\ge 0$ or

• (ii) 

  Stable nodes if and only if $r_H{r}_W>mM$ and ${(r_H-r_W)}^2+4mM\ge 0$.

It is interesting to note that, while the stability of the equilibrium points (or steady state solutions) depends on the emotional inertia ($r_H$ and $r_W$) and the influence function (m and M), their existence depends also on the values of the individual constants $a_H$ and $a_W$, which are proportional to the uninfluenced equilibrium points. For the influence function of the bilinear type, it was discovered that, changing only the slopes m and M, the equilibrium points will move around the quadrants, and their stability will also change. This means that, by changing only the influence functions, the satisfaction of the couple can approach a constant value favorable to both if $(H^{*},W^{*})$ is in quadrant (I), favorable only to the wife if $(H^{*},W^{*})$ is in quadrant (II), unfavorable to both if $(H^{*},W^{*})$ is in quadrant (III), and favorable only to the husband if $(H^{*},W^{*})$ is in quadrant (IV). These steady state solutions can be stable or unstable. If $r_H{r}_W>mM$, the satisfaction of the couple will not change with time; otherwise, it will and can get very far from the values $(H^{*},W^{*})$. If we are given fixed values for the r’s and the a’s, varying the values of the slopes m and M, we can control the position and stability of the equilibrium states. The same is true if the r’s and the slopes are given; then, varying the a’s, stability and position are determined.

In terms of therapeutic interventions, this describes the situation when, without changing each person (the a’s and the r’s), the satisfaction of the couple can be changed by working jointly with both members to modify the values of the interaction (bilinear influence function), or, when working with husband and wife separately, their uninfluenced equilibrium points (the a’s) are changed. On the other hand, the stability of the equilibria can only be altered by changing the relationship between the inertia, $r_W$ and $r_H$, and the slopes $m,M$.

To illustrate the above, we take the values for the variables of the model given in Table 10.1 from Gottman for the validating couple, which are: Inertia for the husband = ${\tilde{r}}_H=0.37$, and inertia for the wife = ${\tilde{r}}_W=0.14$; therefore, $r_H={\tilde{r}}_H-1=-0.63,r_W={\tilde{r}}_W-1=0.86.$ We chose influence functions with $m_{-}=0.31,m_{+}=0.21,$$M_{-}=0.28,$$M_{+}=0.15$. Then, the Jacobians have the form:

$\left(\begin{array}{cc}-0.63\hfill & M\hfill \\ m\hfill & -0.86\hfill \end{array}\right),$

where the values of the slopes $m,M$ depend on the quadrant. For the first quadrant: $m=m_{+}=0.21,M=M_{+}=0.15$. The equilibrium points are

$H^{*}=\frac{0.15a_W+0.86a_H}{0.5103},W^{*}=\frac{0.21a_H+0.63a_W}{0.5103},$

and

$\Delta =0.0529+4\left(0.0315\right)=0.1789>0,r_H{r}_W=0.5418>0.0315=m_{+}{M}_{+};$

then, from Theorem 1, the equilibrium points, if they exist, are stable nodes.

Conditions for the existence of equilibrium points in the first quadrant are $H^{*}>0$ and $W^{*}>0$, which mean:

$0.15a_W+0.86a_H>0\mathrm{and}0.21a_H+0.63a_W>0.$

Similar analyses can be carried out for the other quadrants, which will end up in inequalities between the a’s. In this particular example, for all quadrants, $r_H{r}_W>mM$ and $\Delta >0$; therefore, Theorem 1 implies that all the equilibrium points are stable nodes.

Note that, for general validating personality types, $r_H$ and $r_W$ are both negative, so one of the eigenvalues is negative, but the second one may be positive for large enough values of the radicand. This is consistent with Theorem 1.

For the case of influence functions given by piece-wise constants (11) and (12), the number and values of the equilibrium points depend on the values of the constants, but all equilibrium points have the same local stability.

Theorem 2.

For validating personalities, all the equilibrium points for influence functions given by piece-wise constants (11) and (12) are stable.

The proof is straightforward since the influence functions only add a constant to each equation, and the partial derivatives with respect to H and W are zero, except at the jumps.

Therefore, the Jacobian matrix is

$\left(\begin{array}{cc}{r}_W\hfill & 0\hfill \\ 0\hfill & r_H\hfill \end{array}\right),$

with eigenvalues ${\lambda }_1=r_W,{\lambda }_2=r_H$, and both are negative for validating personalities.

For the other types of influence functions studied, due to the complexity, we will only give the equilibrium points and their stability for fixed values of the parameters. This will be done in the next section.

3. Numerical Simulations

Numerical calculations were done to both validate the theoretical results of the previous section and to find the equilibria and stability for the cases when the analytical solutions are very complicated to find. Particular solutions are also calculated. As mentioned before, we only consider validating personality types since our main objective is to show that the influence functions used have a determining effect on the number of equilibrium points and on their stability and that the same is true of the stochasticity. A similar study can be done for other personality types. The equilibrium points and their stability were found using Mathematica. The simulations of (4) are done using the numerical routine lsode [11], as implemented in Octave. The stochastic system (17) is discretized using Milstein method [9,12].

We are not doing a complete statistical study of the effects of the noise but only including a few representative realizations which show that external influences sometimes do not change very much the state of the system, but sometimes noise can switch the system to a neighborhood of another equilibrium point, with different values for the satisfaction of the couple.

3.1. Single Individual with Validating Personality

First, consider a single individual with validating personality with no external noise. The values of the parameters are chosen to give the same single equilibrium point given by Reference [7]. These are $a_H=0.2394,r_H=-0.63$, which give the uninfluenced equilibrium value of 0.38 for the husband, and $a_W=0.4472,r_W=-0.86$ giving an equilibrium of 0.52 for the wife. Since the eigenvalues are $r_H$ and $r_W$, both uninfluenced equilibrium points are stable. Figure 5 is a simulation of the deterministic and stochastic models for the husband without influence from the wife. For the wife, we get a similar plot. The stochastic simulations were done with $s_1=s_2=0.2$.

3.2. Validating Wife and Husband with Bilinear Influence Function

Next, we consider the influence of the spouse. The simplest case is that of bilinear influence functions (9) and (10). The assumption is that, the greater the happiness, the more influence it has on the spouse. For positive values of the happiness, we used the values given by Reference [7]: $a_H=0.2394,a_W=0.4472,r_H=-0.63,r_w=-0.86$, $m_{+}=0.1922;$$m_{-}=1.7,M_{+}=0.1075,M_{-}=1.9$. The values of the slopes for negative happiness were just chosen to be larger than the ones corresponding to positive happiness, as stated by Reference [7]. The equilibrium states are $(H=0.4873$, $W=0.6289)$ and $(H=-0.3927$, $W=-0.2562)$. The corresponding eigenvalues are ${\lambda }_1=-0.9291,{\lambda }_2=-0.5609$, and ${\lambda }_1=-2.5459,$${\lambda }_2=1.00559$, so the positive equilibrium is stable, and the negative unstable. Figure 6 has a simulation for both the deterministic model and a realization for the stochastic one (top) and a plot of the nullclines showing the points of intersection (bottom). The nullclines are obtained from Equation (4). Setting the first equation equal to zero, we obtain

$a_H+r_{H}H+I_{WH}\left(W\right)=0;$

therefore, one nullcline is given by

$H=H\left(W\right)=\frac{1}{r_H}(-I_{WH}\left(W\right)-a_H).$

Similarly, the other nullcline is given by

$W=W\left(H\right)=\frac{1}{r_W}(-I_{HW}\left(H\right)-a_W).$

Observe that the nullclines are rescaled and displaced versions of the influence functions.

3.3. Piece-Wise Influence Function

We will use the same values of $a_H,a_W,r_H$, and $r_W$ in all our examples. In addition, in all the following examples, we will only consider that $I_{HW}(·)=I_{WH}(·)$ and will call both of them just $I\left(P\right)$ for simplicity, with $P=H$ or $P=W$.

For piece-wise constant influence functions (11) and (12), we studied three different scenarios. Scenario 1 is the case studied by Reference [7], in which the functions are non-decreasing and, as happiness increases, the influence function increases. The values used are $I_{-2}=-2$, $I_{-1}=-1$, $I_1=0.8$, $I_2=1$, $P_{-2}=-1.5$, $P_{-1}=P_1=0$, $P_2=1.5.$ For these values, there are three equilibrium points: $(H=-2.7946,W=-1.8056)$, ($H=-1.2073$, $W=-0.642791)$, and $(H=1.9673,W=1.6828)$. The eigenvalues for the three are ${\lambda }_1=-0.86,{\lambda }_2=-0.63$, so all the equilibrium points are stable. Figure 2 has a simulation for both the deterministic model and a realization for the stochastic one (top) and a plot of the nullclines showing the points of intersection (bottom).

Scenario 2 is similar to scenario 1, except that we consider that there is a region around zero satisfaction for which there is no influence on the spouse. The values of the parameters used are $I_{-2}=-2$, $I_{-1}=-1$, $I_1=0.8$, $I_2=1$, $P_{-2}=-1.5$, $P_{-1}=0.5$, $P_1=0.6$, $P_2=1.5.$ There are now four equilibrium points, $(H=-2.7946,W=-1.8056),(H=-1.2073,$$W=-0.6428),(H=0.38,W=0.52)$, and $(H=1.9673,W=1.6828)$. The eigenvalues are the same as for scenario 1, and all equilibrium points are still stable. Figure 3 has a simulation for both the deterministic model and a realization for the stochastic one (left) and a plot of the nullclines showing the points of intersection (right).

In scenario 3, we take scenario 2, but, instead of using non-decreasing influence functions, we assume that there is a saturation. For large absolute values of the satisfaction, the person pays less attention to the spouse, and their own satisfaction changes slower or even stops changing. Here, we assume that it stops changing. The values of the parameters are $I_{-2}=0$, $I_{-1}=-1$, $I_1=0.8$, $I_2=0$, $P_{-2}=-1.5$, $P_{-1}=0.5$, $P_1=0.6$, $P_2=1.5.$ There are now only two equilibrium points, $(H=-1.2073,W=-0.6428)$ and $(H=0.38,W=0.52)$. The eigenvalues are the same as for scenario 1, and all equilibrium points are still stable. Figure 4 has a simulation for both the deterministic model and a realization for the stochastic one (top) and a plot of the nullclines showing the points of intersection (bottom).

3.4. Piece-Wise Linear Influence Function

The next influence function we consider is a continuous function formed by piece-wise linear functions. For small positive values of the satisfaction, the influence grows linearly up to a point where its starts to decrease linearly and then goes to zero, and similarly for negative values of the satisfaction. The parameters used in the simulations are $P_1=0.5$, $P_2=1.5$, $P_{-1}=-0.5$; $P_{-2}=-1.5$, $I_1=0.7$, $I_2=-1.5.$ There are five equilibrium points, $(H=-1.1405,W=-0.4760)$, $(H=-0.5708,W=-1.1006)$, $(H=-1.0319,$$W=-0.2965)$, $(H=-0.1829,W=-0.1182)$, and $(H=1.1746,W=0.7848)$. For these equilibrium points, the first, third, and fifth are stable, and the other two unstable. Figure 7 has a simulation for both the deterministic model and a realization for the stochastic one with initial values $H_0=0,W_0=0$ (top) and a plot of the nullclines showing the points of intersection (bottom).

Figure 8 shows a simulation for both the deterministic model and a realization for the stochastic one with initial values $H_0=-0.47,W_0=-1.14$. Note that the stochastic solution moves around a different equilibrium point than the deterministic solution. This illustrates that, in certain situations, it is possible for external factors to change the satisfaction of the couple.

The last influence function used in the examples is a piece-wise cubic, continuous and differentiable everywhere. For small positive values of the satisfaction, the influence grows as a cubic, up to a point where its starts to also decrease as a cubic and then goes to zero, and similarly for negative values of the satisfaction. The parameters used in the simulations are $P_1=0.5$, $P_2=1.5$, $P_{-1}=-0.5$;$P_{-2}=-1.5$, $I_1=0.7$, $I_2=-1.5.$ There are seven equilibrium points, $(H=-0.4470,W=-0.8574),(H=-0.6063,$$W=-0.8073),$$(H=-0.7506,$$W=-0.3504),(H=-0.3452,W=-0.2828),(H=1.4651,W=0.5219),(H=0.9063,$$W=0.9179)$, and $(H=0.4787,W=1.2738)$. The first, third, fifth, and seventh are stable, and the other three unstable. Figure 9 has a simulation for both the deterministic model and a realization for the stochastic one with initial values $H_0=0,W_0=0$ (left) and a plot of the nullclines showing the points of intersection (right).

Figure 10 shows a simulation for both the deterministic model and a realization for the stochastic one with initial values $H_0=-1,W_0=-1$. Note that the deterministic solution oscillates between two equilibrium points and that the stochastic solution moves around a different equilibrium point. This illustrates that, in certain situations, it is possible to have periodic deterministic solutions and also for external factors to change the satisfaction of the couple.

4. Conclusions

The type of influence functions used in modeling the satisfaction of a couple has a critical influence on the type of solutions obtained. Even the two types of influence functions considered in Reference [7] produce different number of equilibrium points. This suggests that there is a need for determining those functions more accurately. In addition, external factors, such as economic situation, weather, family, and also more drastic ones, such as illness or a death, unemployment, or getting a better job, can change the satisfaction of a couple. We mostly did concrete examples for a type of personalty, but similar results should happen for different influence functions, parameter values, and personality types.

The analysis of the deterministic model led us to unravel the role of the parameters of the models. While the stability of the equilibrium points (or steady state solutions) depends on the emotional inertia $(r_H$ and $r_W$) and the influence function (m and M), their existence and position depends also on the values of the individual constants $a_H$ and $a_W$. In therapeutic terms, this means that the position of the steady states can be changed by giving therapy separately to husband and wife, but the stability can only be changed by modifying the interaction of the couple.

The study of the stochastic model showed that some perturbations have very little effect on the satisfaction of the couple, but some can move the couple from one steady state to another or to an oscillatory solution. In other words, the introduction of random external effects can sometimes change the satisfaction of the couple, even when the personalities and type of interaction remain unchanged. In our numerical experiments, roughly 20% of cases led to a qualitative change. This is related to the relative sizes of the white noise and the basins of attraction of the different equilibrium points.

For ease of reference, we include, in Table 1, a summary of the data for the numerical results.

Author Contributions

Conceptualization, B.C.-C., C.E.G.-H. and M.d.C.J.; methodology, B.C.-C., C.E.G.-H. and M.d.C.J.; software, B.C.-C., C.E.G.-H. and M.d.C.J.; validation, B.C.-C., C.E.G.-H. and M.d.C.J.; formal analysis, B.C.-C., C.E.G.-H. and M.d.C.J.; investigation, B.C.-C., C.E.G.-H. and M.d.C.J.; resources, B.C.-C., C.E.G.-H. and M.d.C.J.; data curation, B.C.-C., C.E.G.-H. and M.d.C.J.; writing—original draft preparation, B.C.-C., C.E.G.-H. and M.d.C.J.; writing—review and editing, B.C.-C., C.E.G.-H. and M.d.C.J.; visualization, B.C.-C., C.E.G.-H. and M.d.C.J.; supervision, B.C.-C., C.E.G.-H. and M.d.C.J.; project administration, B.C.-C., C.E.G.-H. and M.d.C.J. All authors have read and agreed to the published version of the manuscript.

Funding

This research received no external funding.

Acknowledgments

We thank Ramiro Chavez and Ana Perez for their computational support and the reviewers for their very pertinent comments.

Conflicts of Interest

The authors declare no conflict of interest.

Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Figures and Table

Figure 1. Bilinear influence functions IHW(H) (top) and IWH(W) (bottom).

Figure 2. Piece-wise constant influence functions, scenario 1. Plot of H and W vs. time for the deterministic case and a stochastic realization (top) and a plot of the nullclines and the equilibrium points (bottom). Observe that these nullclines are discontinuous; therefore, for the H nullcline only, the horizontal lines count, and, for the W nullcline, only the vertical parts count.

Figure 3. Piece-wise constant influence functions, scenario 2. Plot of H and W vs. time for the deterministic case and a stochastic realization (top) and a plot of the nullclines and the equilibrium points (bottom).

Figure 3. Piece-wise constant influence functions, scenario 2. Plot of H and W vs. time for the deterministic case and a stochastic realization (top) and a plot of the nullclines and the equilibrium points (bottom).

Figure 4. Piece-wise constant influence functions, scenario 2. Plot of H and W vs. time for the deterministic case and a stochastic realization (top) and a plot of the nullclines and the equilibrium points (bottom). In this example, the external white noise significantly changes the state of satisfaction of the couple.

Figure 5. A numerical simulation of the deterministic case and a stochastic realization for a validating husband with no influence from the wife.

Figure 6. Bilinear influence functions. Plot of H and W vs. time for the deterministic case and a stochastic realization (top) and a plot of the nullclines and the equilibrium points (bottom). In this example, the external noise only makes the satisfaction oscillate around the equilibrium state of the deterministic case.

Figure 7. Piece-wise linear influence functions. Plot of H and W vs. time for the deterministic case and a stochastic realization with initial values H0=0,W0=0 (top) and a plot of the nullclines and the equilibrium points (bottom).

Figure 8. Piecewise linear influence functions. Plot of H and W vs. time for the deterministic case and a stochastic realization with initial values H0=−0.47, W0=−1.14.

Figure 9. Piece-wise cubic influence functions. Plot of H and W vs. time for the deterministic case and a stochastic realization with initial values H0=0,W0=0 (left) and a plot of the nullclines and the equilibrium points (right).

Figure 10. Piece-wise cubic influence functions. Plot of H and W vs. time for the deterministic case and a stochastic realization with initial values H0=−1,W0=−1.

Table 1

Summary of parameter values used. The same values for $a_H,a_W,r_H$, and $r_W$ are used in all our examples.

Figure Number	Type of Influence Function	Parameters
2	Zero influence function	$a_H=0.2394$, $r_H=-0.63$.
3	Bilinear	$a_H=0.2394$, $a_W=0.4472$, $r_H=-0.63$, $r_w=-0.86$, $m_{+}=0.1922$; $m_{-}=1.7$, $M_{+}=0.1075$, $M_{-}=1.9$.
4	Piece-wise constant	$I_{-2}=-2$, $I_{-1}=-1$, $I_1=0.8$, $I_2=1$, $P_{-2}=-1.5$, $P_{-1}=P_1=0$, $P_2=1.5.$
5	Piece-wise constant	$I_{-2}=-2$, $I_{-1}=-1$, $I_1=0.8$, $I_2=1$, $P_{-2}=-1.5$, $P_{-1}=0.5$, $P_1=0.6$, $P_2=1.5.$
6	Piece-wise constant with saturation	$I_{-2}=0$, $I_{-1}=-1$, $I_1=0.8$, $I_2=0$, $P_{-2}=-1.5$, $P_{-1}=0.5$, $P_1=0.6$, $P_2=1.5.$
7 and 8	Piece-wise linear	$P_1=0.5$, $P_2=1.5$, $P_{-1}=-0.5$; $P_{-2}=-1.5$, $I_1=0.7$, $I_2=-1.5.$
9 and 10	Piece-wise cubic	$P_1=0.5$, $P_2=1.5$, $P_{-1}=-0.5$; $P_{-2}=-1.5$, $I_1=0.7$, $I_2=-1.5.$