Furthermore, the iteration on vertices might provide a new tool for an
empirical vertex density analysis in the case where the training loss can be
written as a ReLU feed-forward neural network. For example as discussed
previously in Example 22, the loss decay behaviour there seems to behave
exponential in the number of visited vertices such that this training loss
curve might contain information on the quality of a local minimum converged
to. In addition, local vertex density could be estimated by our DRLSimplex
algorithm using the distance of vertex positions which in turn could serve as
an input to tune or improve the step size control in other gradient descent
training algorithms.

##  5 Summary

We have introduced the class of deep ReLU programming problems as non-convex
optimization problems of finding a local minimum of a real-valued feed-forward
neural network $f$ with ReLU activation functions. These functions are piece-
wise affine on domains which are given by sets of linear inequalities. In
particular, if restricted to one such domain, the resulting minimization
problem has the form of a linear programming problem where the restricted
domain takes the role of the feasible region. Hence, deep ReLU programming can
be seen as a generalization of linear programming with the difference that we
allow multiple feasible regions, each with its own affine objective function
such that the overall function is continuous.

We developed notation and theory suitable for the analysis of these functions
in Section 2. For an input value $x$ we defined activation patterns $s$ to be
compatible if replacing the ReLU activation functions by hard-coded activities
specified in $s$ preserves the output and all intermediate layer values of
$f$, i.e. if the network $f$ can be described by fixed affine layer maps
induced by $s$. We call quantities _objective_ when they are defined based on
the true network $f$ and _subjective_ when they are defined based on this
fixed affine layer network induced by some activation pattern $s$. Quantities
that have an objective and meaningful subjective versions are _stable_ , if
the objective version coincides with the subjective versions for compatible
$s$. We have introduced several local network quantities such as the critical
kernel and the critical indices and proved stability. We further formulated
the concept of _regularity_ which is a requirement in our stability result for
the critical indices. For a random continuous independent distribution on its
parameters the resulting network is almost surely regular which shows that
regularity is not a strong requirement. Our stability results lay the
foundation for algorithms that use the activation of the neurons as a state
parameters because for they identify subjective quantities that do not need to
be updated when the state parameter $s$ is changed from one compatible
activation pattern to another, thus a change of the “feasible region”
identified by $s$ without the need to update or recompute such quantities. We
then formalized the notion of a vertex, its corresponding axes that separate
affine regions of $f$ and formulated a necessary and sufficient condition for
$x$ being a local minimum of $f$.

In addition we provided algorithmic subroutines suitable for the construction
of iterative optimization algorithms for deep ReLU programming problems in
Section 3. In Section 4.1 we combined these subroutines with our theoretical
concepts and results to formalize our DRLSimplex algorithm, a generalization
of the standard simplex algorithm to deep ReLU programming that can be used to
iteratively optimize regular feed-forward ReLU neural networks. In its first
phase it first finds a vertex and then iterates on vertices while changing the
affine region and keeping track of a compatible activation pattern and a set
of critical indices that define the region separating axes to proceed with.
The automatic step size selection exploits the piece-wise constant gradient
and ensures that the objective function cannot increase across iterations. Our
algorithm stops based on the local minimum condition for regular networks we
formulated previously. For an $L$ layer network of widths $n_{0},\ldots,n_{L}$
we showed that every iteration is of computational complexity
$\mathcal{O}(\sum_{i=1}^{L}n_{i}n_{i-1}+n_{0}^{2})$, where the first and the
second terms reflect the number of parameters and simplex tableau size, such
that this complexity order is optimal for an extension of the simplex
algorithm to ReLU neural networks. We further proposed modifications for
better numeric stability, other activation functions, alternatives to vertex
iterations and generalizations for mixtures of feed forward ReLU nets and
quadratic forms as objective functions.

For our DRLSimplex algorithm we provided a proof-of-concept implementation in
the Julia programming language in Section 4.2 which is demonstrated in Section
4.3 at optimization problems that can be rewritten in the form of a deep ReLU
programming problem. We showed that our algorithm can be applied in LASSO
penalized quantile regression, censored least absolute deviation and even for
training of the first layer of a deep ReLU network using L1 loss. For the
first and second applications in this list, there exist specialized linear
programming solver algorithms that can be seen as special cases of our
DRLSimplex algorithm, for the third example one usually resorts to standard
neural network training algorithms not specifically designed to take advantage
of the piece-wise affine structure of the ReLU network to be trained. Here our
algorithm provides a simplex-like iteration on vertices of the loss function
such that the decreases in every step. In future work it will be interesting
to apply deep linear programming for training of not only the first layer of
ReLU neural networks but also other layers iteratively. This would yield a
training procedure which monotonically decreases the L1 training loss across
iterations, a major advantage compared to other optimization methods. We also
demonstrated how our DRLSimplex method can be generalized for objective
functions that are a mixture of a feed-forward ReLU neural network and a
positive definite form: Our implementation is able to successfully perform
LASSO optimization.

We finally discussed our DRLSimplex algorithm and compared it to linear
programming and gradient descent variants in Section 4.4. The overall
computational complexity for an iteration in our algorithm of
$\mathcal{O}(n_{0}^{2}+\sum_{i=1}^{L}n_{i}n_{i-1})$ allows it to be considered
as a new competitive optimization tool. Especially for loss functions that can
be rewritten in the form of a feed-forward ReLU neural network, in a
combination other gradient-descent based algorithms could benefit from the
strictly decreasing loss guaranteed by our algorithm’s automatic step size
control or the fact that exact vertex positions are returned. In particular
the latter property could be used to empirically analyze the vertex density of
such loss functions which in turn could provide new insight in their structure
or allow for improved step-size control of gradient-descent algorithms based
on local vertex density.

We hope that our theoretical framework, our algorithmic building blocks and
our DRLSimplex algorithm open the door for the application of linear
programming and active-set related techniques to feed-forward ReLU neural
networks. This could help to analyze and better understand the vertex
structure of loss functions related to ReLU neural networks and might
eventually lead to improved neural network training algorithms.

##  Appendix A Proofs and auxiliary results

In the sequel, we use the notation introduced in Sections 1 and 2.

###  A.1 Compatibility

###### Proof of Theorem 2.

Fix $s\in\mathcal{S}^{C}(x)$. For $(l,j)\in\mathcal{I}$, if $s_{lj}=0$ then by
equation (14) $\textnormal{sign}(A(x)_{lj})=H(x)_{lj}\in\left\\{0,-1\right\\}$
and therefore $\textnormal{ReLU}(A(x)_{lj})=0$. Similarly $s_{lj}=1$ implies
$\textnormal{ReLU}(A(x)_{lj})=A(x)_{lj}$. Hence
$g^{(1)}(x)=\tilde{g}^{(1)}_{s_{1}}(x)$ and for
$l\in\left\\{2,\ldots,L\right\\}$,

$g^{(l)}\circ\cdots\circ
g^{(1)}(x)=\textnormal{ReLU.}(A(x)_{l})=\textnormal{diag}(s_{l})A(x)_{l}=\tilde{g}^{(l)}_{s_{l}}\circ
g^{(l-1)}\circ\cdots\circ g^{(1)}(x).$

A recursive application of the above formula yields $A(x)=\tilde{A}_{s}(x)$.

Now fix $s\in\tilde{\mathcal{S}}^{C}(x)$. By equation (15) $s$ is compatible
with $\tilde{H}_{s}(x)$. In particular for $(l,j)\in\mathcal{I}$, $s_{lj}=0$
implies
$\textnormal{sign}(\tilde{A}_{s}(x)_{lj})=\tilde{H}_{s}(x)_{lj}\in\left\\{-1,0\right\\}$
and $\textnormal{ReLU}(\tilde{A}_{s}(x)_{lj})=0$. Similarly $s_{lj}=1$ implies
$\textnormal{ReLU}(\tilde{A}_{s}(x)_{lj})=\tilde{A}_{s}(x)_{lj}$. This implies

$\forall
l\in\left\\{1,\ldots,L\right\\}\quad\textnormal{ReLU.}(\tilde{A}_{s}(x)_{l})=\textnormal{diag}(s_{l})\tilde{A}_{s}(x)_{l}$
(52)

We now prove $A(x)=\tilde{A}_{s}(x)$ by induction. Note that
$A(x)_{1}=W_{1}x+b_{1}=\tilde{A}_{s}(x)_{1}$. For the induction step assume
$A(x)_{l}=\tilde{A}_{s}(x)_{l}$ for some $l\in\left\\{1,\ldots,L-1\right\\}$.
Then by equation (52)

$g^{(l)}\circ\cdots\circ
g^{(1)}(x)=\textnormal{ReLU.}(A(x)_{l})=\textnormal{ReLU.}(\tilde{A}_{s}(x)_{l})=\textnormal{diag}(s_{l})(\tilde{A}_{s}(x)_{l})=\tilde{g}_{s_{l}}^{(l)}\circ\cdots\circ\tilde{g}_{s_{1}}^{(1)}(x).$

This implies that $A(x)_{l+1}=W_{l+1}g^{(l)}\circ\cdots\circ
g^{(1)}(x)+b_{l+1}=W_{l+1}\tilde{g}_{s_{l}}^{(l)}\circ\cdots\circ\tilde{g}_{s_{1}}^{(1)}+b_{l+1}=\tilde{A}_{s}(x)_{l+1}$.
An induction on $l$ yields $A(x)=\tilde{A}_{s}(x)$.

Therefore we have shown that
$s\in\mathcal{S}^{C}(x)\cup\tilde{\mathcal{S}}^{C}(x)$ implies
$A(x)=\tilde{A}_{s}(x)$, hence $H(x)=\tilde{H}_{s}(x)$ and therefore
$\mathcal{S}^{C}(x)=\tilde{\mathcal{S}}^{C}(x)$. ∎

###### Proof of Corollary 3.

This follows from Theorem 2 because
$\tilde{f}_{s}(x)=W_{L+1}\textnormal{ReLU}.(\tilde{A}_{s}(x)_{L})+b_{L+1}$ and
$f(x)=W_{L+1}\textnormal{ReLU}(A(x)_{L})+b_{L+1}$. ∎

###### Proof of Corollary 4.

Let $y\in\mathbb{R}^{n_{0}}$ with $H(y)=H(x)$. Then $H(y)$ and $H(x)$ are
compatible with $s$ such that Theorem 2 implies that
$\tilde{H}_{s}(y)=\textnormal{sign.}(\tilde{A}_{s}(y))=\textnormal{sign.}(A(y))=H(y)$
and similarly $\tilde{H}_{s}(x)=H(x)$. In particular
$\tilde{H}_{s}(y)=\tilde{H}_{s}(x)$. For the other direction assume
$y\in\mathbb{R}^{n_{0}}$ with $\tilde{H}_{s}(y)=\tilde{H}_{s}(x)$. By
assumption and again Theorem 2 it holds that $s\in\tilde{\mathcal{S}}^{C}(x)$,
i.e. $s$ is compatible with $\tilde{H}_{s}(x)$ and hence also with
$\tilde{H}_{s}(y)$ such that $s\in\tilde{\mathcal{S}}^{C}(y)$. Like above, it
follows from the same theorem that
$H(y)=\tilde{H}_{s}(y)=\tilde{H}_{s}(x)=H(x)$. ∎

###  A.2 Local Behaviour

######  Lemma 23.

For all $x\in\mathbb{R}^{n_{0}}$, there exists $\varepsilon>0$ such that for
all $y\in B_{\varepsilon}(x)$ it holds that
$\mathcal{S}^{C}(y)\subset\mathcal{S}^{C}(x)$.

###### Proof.

Let $\varepsilon>0$ such that $H(x)_{lj}=H(y)_{lj}$ for $y\in
B_{\varepsilon}(x)$ and $(l,j)\in\mathcal{I}\setminus C(x)$. By Lemma 5
$H(x)_{lj}=0$ for $(l,j)\in C(x)$ such that for all $y\in B_{\varepsilon}(x)$
any $s\in\mathcal{S}$ compatible with $H(y)$ is also compatible with $H(x)$. ∎

This means that for every input $x\in\mathbb{R}^{n_{0}}$ there exists a
neighbourhood such that for all $y$ in this neighbourhood all $s$ compatible
with $H(y)$ are also compatible with $H(x)$.

###### Proof of Lemma 9.

If $v^{*}\in\tilde{\textnormal{Ker}}{}^{C}_{s}(x)$, then there exists
$\varepsilon^{*}>0$ such that for all $(l,j)\in\tilde{C}_{s}(x)$ and
$t\in(-\varepsilon,\varepsilon)$
$\tilde{H}_{s}(x+tv^{*})_{lj}=\tilde{H}_{s}(x)_{lj}=0$. In particular for all
$t\in(-\varepsilon,\varepsilon)$ it holds that
$\tilde{A}_{s}(x+tv^{*})_{lj}=0$. By affinity of
$z\mapsto\tilde{A}_{s}(z)_{lj}$, it holds that
$\tilde{A}_{s}(x+v^{*})_{lj}=0$. For the other direction assume that
$v^{*}\in\mathbb{R}^{n_{0}}$ with $\tilde{A}_{s}(x+v^{*})_{lj}=0$ for
$(l,j)\in\tilde{C}_{s}(x)$. Let $\varepsilon^{*}>0$ such that for all $v\in
B_{\varepsilon^{*}}(0)$ and $(l,j)\mathcal{I}\setminus\tilde{C}_{s}(x)$ it
holds that $\tilde{H}_{s}(x+v)_{lj}=\tilde{H}_{s}(x)_{lj}$. This is possible
by equation (17). Then for all $t\in\mathbb{R}^{n_{0}}$ such that $tv^{*}\in
B_{\varepsilon^{*}}(0)$ it holds that
$\tilde{H}_{s}(x+tv^{*})=\tilde{H}_{s}(x)$ since also for indices
$(l,j)\in\tilde{C}_{s}(x)$ it holds that
$\tilde{H}_{s}(x+tv)=0=\tilde{H}_{s}(x)$ by Lemma 5. ∎

###### Proof of Lemma 10.

By Corollary 4 $\big{\\{}v\in\mathbb{R}^{n_{0}}\mid
H(x+v)=H(x)\big{\\}}=\big{\\{}v\in\mathbb{R}^{n_{0}}\mid\tilde{H}_{s}(x+v)=\tilde{H}_{s}(x)\big{\\}}$.
Let $\varepsilon^{*}>0$ small enough such that $H_{s}(x+v)_{ij}$ is constant
for indices $(i,j)\in\mathcal{I}\setminus\tilde{C}_{s}(x)$ and $v\in
B_{\varepsilon^{*}}(0)$. This is possible by equation (17). For all other
indices $(i,j)\in\tilde{C}_{s}(x)$, the subjective hyperplane pattern
satisfies $\tilde{H}_{s}(x)_{ij}=0$ by Lemma 5. It follows that
$B_{\varepsilon^{*}}(x)\cap\big{\\{}v\in\mathbb{R}^{n_{0}}\mid\tilde{H}_{s}(x+v)=\tilde{H}_{s}(x)\big{\\}}=B_{\varepsilon^{*}}(x)\bigcap_{(l,j)\in\tilde{C}_{s}(x)}\left\\{v\in\mathbb{R}^{n_{0}}\mid\tilde{A}_{s}(x+v)_{lj}=0\right\\}$.
Lemmas 9 and 8 complete the proof. ∎

###  A.3 Regularity

####  A.3.1 Basic results

######  Lemma 24.

Assume that $x\in\mathbb{R}^{n_{0}}$, $(l^{*},j^{*})\in\mathcal{I}$ and
$s,s^{\prime}\in\mathcal{S}^{C}(x)$ with $s_{lj}=s^{\prime}_{lj}$ for
$(l,j)\in\mathcal{I}\setminus\left\\{(l^{*},j^{*})\right\\}$ and
$s_{l^{*}j^{*}}\neq s^{\prime}_{l^{*}j^{*}}$. If $v_{s,\tilde{l},\tilde{j}}=0$
and $v_{s^{\prime},\tilde{l},\tilde{j}}\neq 0$ for some
$(\tilde{l},\tilde{j})\in\mathcal{I}$ then
$v_{s^{\prime},\tilde{l},\tilde{j}}=\alpha v_{s,l^{*},j^{*}}$ for some
$\alpha\in\mathbb{R}$.

###### Proof.

By equation (24) $\tilde{l}>l^{*}$ and

$v_{s^{\prime},\tilde{l},\tilde{j}}=v_{s^{\prime},\tilde{l},\tilde{j}}-v_{s,\tilde{l},\tilde{j}}=\left(W_{\tilde{l}}\textnormal{diag}(s_{\tilde{l}-1})W_{l-1}\cdots\textnormal{diag}(s_{l^{*}+1})W_{l^{*}+1}\right)_{\tilde{j}j^{*}}v_{s,l^{*},j^{*}}$

for $s^{\prime}_{j^{*}l^{*}}=1$. Similarly, for $s^{\prime}_{l^{*}j^{*}}=0$ we
obtain the negative right-hand side. ∎

###### Proof of Theorem 14.

For two compatible activation patterns $s,s^{\prime}\in\mathcal{S}^{C}(x)$
that differ only at one position $(l^{*},j^{*})\in\mathcal{I}$, it holds that
$(l^{*},j^{*})\in\tilde{C}_{s}(x)\cap\tilde{C}_{s^{\prime}}(x)$ by equations
(11) and (17). If $\tilde{C}_{s}(x)\neq\tilde{C}_{s^{\prime}}(x)$ then there
exists
$(l^{\prime},j^{\prime})\in\tilde{C}_{s}(x)\cup\tilde{C}_{s^{\prime}}(x)$ with
both, $v_{s,l^{*},j^{*}}=0$ and $v_{s^{\prime},l^{*},j^{*}}\neq 0$ or both,
$v_{s^{\prime},l^{*},j^{*}}=0$ and $v_{s,l^{*},j^{*}}\neq 0$. But by Lemma 24
this violates the regularity assumption such that
$\tilde{C}_{s}(x)=\tilde{C}_{s^{\prime}}(x)$.

For two arbitrary compatible activation patterns
$s,s^{\prime}\in\mathcal{S}^{C}(x)$ that differ at $k>1$ positions we can find
a transformation path
$s=:s^{(1)},\ldots,s^{(k+1)}:=s^{\prime}\in\mathcal{S}^{C}(x)$ such that two
subsequent activation patterns $s^{(i)},s^{(i+1)}$,
$i\in\left\\{1,\ldots,k\right\\}$ differ only at one position
$(l_{i},j_{i})\in\mathcal{I}$. Now the above result implies
$\tilde{C}_{s^{(1)}}(x)=\ldots=\tilde{C}_{s^{(k)}}(x)$, in particular
$\tilde{C}_{s}(x)=\tilde{C}_{s^{\prime}}(x)$.

We now have shown that there exists $C\subset\mathcal{I}$ with
$\tilde{C}_{s}(x)=C^{*}$ for all $s\in\mathcal{S}^{C}(x)$. Now equations (16)
and $\eqref{eq:subjectiveC}$ imply that this set $C^{*}$ of indices must be
exactly $C(x)$. ∎

####  A.3.2 Probabilistic results on regular vertices

Below we use the notation and concepts of Section 2. For $s\in\mathcal{S}$,
$(l,j)\in\mathcal{I}$ let

$H^{s}_{lj}:=\left\\{x\in\mathbb{R}^{n_{0}}\middle|\;\langle
x,\tilde{v}_{s,l,j}\rangle+\left(\sum_{i=1}^{l-1}W_{l}\textnormal{diag}(s_{l-1})\cdots
W_{i+1}\textnormal{diag}(s_{i})b_{i}\right)_{j}+b_{lj}=0\right\\}.$ (53)

Note that
$H^{s}_{lj}=\big{\\{}x\in\mathbb{R}^{n_{0}}\mid\;\tilde{H}_{s}(x)_{lj}=0\big{\\}}$.

######  Lemma 25.

Assume the $\sum_{i=2}^{L+1}n_{l}n_{l-1}$ weight parameters of the matrices
$W_{1},\ldots,W_{L}$ and the $\sum_{i=2}^{L+1}n_{i}$ bias parameters of the
vectors $b_{1},\ldots,b_{L}$ of the neural network $f$ are sampled from a
distribution such that conditionally on the weight parameters, the bias
parameters are independent with a conditional marginal distribution that
assigns probability zero to all finite sets. Then the following holds almost
surely: For all $n\in\mathbb{N}$, $s_{1},\ldots,s_{n}\in\mathcal{S}$,
$(l_{1},j_{1}),\ldots,(l_{n},j_{n})\in\mathcal{I}$ with non-zero normal
vectors $\tilde{v}_{s_{1},l_{1},j_{1}},\ldots,\tilde{v}_{s_{n},l_{n},j_{n}}$
the hyperplanes $H^{s_{1}}_{l_{1},j_{1}},\ldots,H^{s_{n}}_{l_{n},j_{n}}$
defined by equation (53) satisfy either
$\cap_{i=1}^{n}H^{s_{i}}_{l_{i},j_{i}}=\left\\{\right\\}$ or the normal
vectors $\tilde{v}_{s_{1},l_{1},j_{1}},\ldots,\tilde{v}_{s_{n},l_{n},j_{n}}$
are linearly independent.

###### Proof.

We proceed by induction. For $n=1$ the claim is obviously true. Now assume the
claim is true for $n=n^{*}\in\mathbb{N}$ and let
$s_{1},\ldots,s_{n^{*}+1}\in\mathcal{S}$,
$(l_{1},j_{1}),\ldots,(l_{n^{*}+1},j_{n^{*}+1})\in\mathcal{I}$. Conditionally
on the weight parameters, the normal vectors

$\tilde{v}_{s_{1},l_{1},j_{1}},\ldots,\tilde{v}_{s_{n^{*}+1},l_{n^{*}+1},j_{n^{*}+1}}$
(54)

are almost surely constant and by assumption the bias parameters are
independent with a continuous distribution. We distinguish two cases:

  1. 1.

Conditionally on the weight parameters, at least one of the above normal
vectors is equal to zero.

  2. 2.

Conditionally on the weight parameters, all of the above normal vectors are
nonzero. In this case let $k\in\left\\{1,\ldots,n^{*}+1\right\\}$ with
$l_{k}=\max_{i\in\left\\{1,\ldots,n^{*}\right\\}}l_{i}$ and let
$J:=\left\\{1,\ldots,n^{*}+1\right\\}\setminus\\{k\\}$. Now we condition on
the weight parameters and on $b_{l_{i},j_{i}}$, $i\in J$. In this case, the
bias vector $b_{l_{k},j_{k}}$ has a continuous distribution. If
$\cap_{i=1}^{n^{*}+1}H^{s_{i}}_{l_{i},j_{i}}\neq\\{\\}$ then by induction
assumption the vectors $\tilde{v}_{s_{i},l_{i},j_{i}}$, $i\in J$ are linearly
independent. If $\tilde{v}_{s_{i},l_{i},j_{i}}$,
$i\in\left\\{1,\ldots,n^{*}\right\\}$ are linearly dependent,
$\cap_{i=1}^{n^{*}+1}H^{s_{i}}_{l_{i},j_{i}}\neq\\{\\}$ implies that
$\cap_{i\in J}H^{s_{i}}_{l_{i},j_{i}}\subset H^{s_{k}}_{l_{k},j_{k}}$ and this
is only possible when $b_{l_{k},j_{k}}$ attains a specific value
$b^{*}\in\mathbb{R}^{*}$ which almost surely does not happen.

The above discussion shows that with probability zero the normal vectors in
equation (54) are all non-zero, linearly dependent and the hyperplanes satisfy
$\cap_{i=1}^{n^{*}+1}H^{s_{i}}_{l_{i},j_{i}}\neq\\{\\}$. Since there are only
finitely many choices for $s_{1},\ldots,s_{n^{*}+1}\in\mathcal{S}$,
$(l_{1},j_{1}),\ldots,(l_{n^{*}+1},j_{n^{*}+1})\in\mathcal{I}$ the claim holds
almost surely for $n=n^{*}+1$. ∎

###### Proof of Theorem 15.

The claim is just a reformulation of Lemma 25. ∎

###  A.4 Feasible directions

######  Lemma 26.

For $x\in\mathbb{R}^{n_{0}}$ it holds that $s\in\mathcal{S}^{C}(x)$ if and
only if $\tilde{A}_{s}(x)_{lj}(s_{lj}-0.5)\geq 0$ for all
$(l,j)\in\mathcal{I}$.

###### Proof.

This follows from equations (13) and (12). ∎

###### Proof of Lemma 16.

Let $v^{*}\in\mathbb{R}^{n_{0}}$ with
$\tilde{A}_{s}(x+v^{*})_{lj}(s_{lj}-0.5)\geq 0$ for all $\tilde{C}_{s}(x)$.
Choose $\varepsilon^{*}>0$ small enough such that for indices
$(l,j)\in\mathcal{I}\setminus\tilde{C}_{s}(x)$, the subjective hyperplane
pattern is constant, i.e. such that for all $v\in B_{\varepsilon^{*}}(0)$ it
holds that $\tilde{H}_{s}(x+v)_{lj}=\tilde{H}_{s}(x)_{lj}$. This is possible
by equation (17). Now let $\varepsilon>0$ such that $\varepsilon v^{*}\in
B_{\varepsilon^{*}}(0)$. In particular at indices
$(l,j)\in\mathcal{I}\setminus\tilde{C}_{s}(x)$ $\tilde{H}_{s}(x+\varepsilon
v^{*})_{lj}=\tilde{H}_{s}(x)_{lj}$ is compatible with $s$. For the other
indices $(l,j)\in\tilde{C}_{s}(x)$ it holds by assumption and Lemma 12 that

$\tilde{A}_{s}(x+\varepsilon v^{*})_{lj}(s_{lj}-0.5)=\varepsilon\langle
v^{*},\tilde{v}_{s,l,j}\rangle(s_{lj}-0.5)=\varepsilon\tilde{A}_{s}(x+v^{*})_{lj}(s_{lj}-0.5)\geq
0.$

Using Lemma 26 it follows that $\tilde{H}_{s}(x+\varepsilon v^{*})$ is
compatible with $s$. To prove the other direction let $v\in\mathbb{R}^{n_{0}}$
such that for some $(l,j)\in\tilde{C}_{s}(x+v)$ the condition
$A_{s}(v)_{lj}(s_{lj}-0.5)\geq 0$ does _not_ hold. Similar to above it follows
that $H(x+\varepsilon v)$ cannot be compatible with $s$ for any
$\varepsilon>0$. ∎

######  Corollary 27.

For an input $x\in\mathbb{\mathbb{missing}}{R}^{n_{0}}$ and a compatible
activation pattern $s\in\mathcal{S}^{C}(x)$ it holds that all directions
$v\in\mathbb{R}^{n_{0}}$ such that $s\in\mathcal{S}^{C}(x+v)$ satisfy
$v\in\tilde{D}_{s}(x)$.

###### Proof of Corollary 27.

This follows immediately from Lemmas 26 and 16. ∎

The above corollary gives a sufficient criterion for directions being feasible
directions in terms of compatibility.

###### Proof of Theorem 18.

For $r:=s_{lj}=s^{\prime}_{lj}$, let

$\mathcal{W}_{s}:=\left\\{v\in\mathbb{R}^{n_{0}}\mid\forall(l,j^{\prime})\in\mathcal{I}\setminus\left\\{(l,j)\right\\}\;\tilde{H}_{s}(x+v)_{l^{\prime}j^{\prime}}=\tilde{H}_{s}(x)_{l^{\prime}j^{\prime}}\land(r-0.5)\tilde{A}_{s}(x+v)_{lj}>0\right\\}$

and define $\mathcal{W}_{s^{\prime}}$ similarly. Corollary 4 implies
$\mathcal{W}_{s}=\mathcal{W}_{\tilde{s}}$. By equation (17) and Theorem 14 it
is possible to find $\varepsilon_{s},\varepsilon_{s^{\prime}}>0$ such that for
$y_{s}\in B_{\varepsilon_{s}}(x)$, $y_{s^{\prime}}\in
B_{\varepsilon_{s^{\prime}}}(x)$ the hyperplane patterns
$\tilde{H}_{s}(y_{s})_{lj}$, $\tilde{H}_{s^{\prime}}(y_{s^{\prime}})_{lj}$ are
constant at indices $(l,j)\not\in
C(x)=\tilde{C}_{s}(x)=\tilde{C}_{s^{\prime}}(x)$. Since
$B_{\min(\varepsilon_{s},\varepsilon_{s^{\prime}})}(x)$ contains an open ball
around $x$, it follows that

$\displaystyle\;\left\\{v\in\mathbb{R}^{n_{0}}\mid\forall(l,j^{\prime})\in\tilde{C}_{s}(x)\setminus\left\\{(l,j)\right\\}\;\tilde{H}_{s}(x+v)_{l^{\prime}j^{\prime}}=\tilde{H}_{s}(x)_{l^{\prime}j^{\prime}}\land(r-0.5)\tilde{A}_{s}(x+v)_{lj}>0\right\\}=$
$\displaystyle\;\left\\{v\in\mathbb{R}^{n_{0}}\mid\forall(l,j^{\prime})\in\tilde{C}_{s^{\prime}}(x)\setminus\left\\{(l,j)\right\\}\;\tilde{H}_{s^{\prime}}(x+v)_{l^{\prime}j^{\prime}}=\tilde{H}_{s^{\prime}}(x)_{l^{\prime}j^{\prime}}\land(r-0.5)\tilde{A}_{s^{\prime}}(x+v)_{lj}>0\right\\}.$

Now note that by equation (28), $\tilde{a}_{s,x,l,j}$ is an element of the
first of the above sets and by equivalence also of the second. In particular,
for $\varepsilon^{*}>0$ small enough
$x+\varepsilon^{*}\tilde{a}_{s,x,l,j}\in\mathcal{W}_{s}=\mathcal{W}_{\tilde{s}}$
such that
$s,s^{\prime}\in\mathcal{S}^{C}(x+\varepsilon^{*}\tilde{a}_{s,x,l,j})$. Now
Lemma 12 and Theorem 2 imply

$\displaystyle\varepsilon^{*}\langle\tilde{a}_{s,x,l,j},\tilde{v}_{s^{\prime},l^{\prime},j^{\prime}}\rangle$
$\displaystyle=\tilde{A}_{s^{\prime}}(x+\varepsilon^{*}a_{s,x,l,j})_{l^{\prime}j^{\prime}}=A(x+\varepsilon^{*}\tilde{a}_{s,x,l,j})_{l^{\prime}j^{\prime}}=\tilde{A}_{s}(x+\varepsilon^{*}\tilde{a}_{s,x,l,j})_{l^{\prime}j^{\prime}}$
$\displaystyle=\varepsilon^{*}\langle\tilde{a}_{s,x,l,j},\tilde{v}_{s,l^{\prime},j^{\prime}}\rangle=\varepsilon^{*}\delta_{(l,j),(l^{\prime},j^{\prime})}.$

for all critical indices $(l^{\prime},j^{\prime})\in C(x)$. By equation (28)
this means $\tilde{a}_{s,x,l,j}=\tilde{a}_{s^{\prime},x,l,j}$. ∎

###  A.5 Regular vertices and local minima

###### Proof of Lemma 21.

If a strict local minimum $x\in\mathbb{R}^{n_{0}}$ is not a vertex then there
exists a non-zero $v\in\textnormal{Ker}^{C}(x)$. But by equation (18) cannot
be a strict local minimum. For the second claim let $x^{*}$ be a non-strict
local minimum of $f$, i.e. there exists $\varepsilon^{*}>0$ such that
$f(x^{*})=\inf\big{\\{}f(x)\mid\;x\in B_{\varepsilon}(x^{*})\big{\\}}$. If
$\textnormal{df}^{C}(x^{*})>0$ then choose a non-zero
$v\in\textnormal{Ker}^{C}(x^{*})$. The map $t\mapsto f(x^{*}+tv)$ is an affine
function for $t\in\left\\{t\in\mathbb{R}\mid\;H(x^{*}+tc)=H(x^{*})\right\\}$.
By equation (18) there exists $\varepsilon>0$ such that
$(-\varepsilon,\varepsilon)\subset\left\\{t\in\mathbb{R}\mid\;H(x^{*}+tv)=H(x^{*})\right\\}$.
This implies that, the function $\mathbb{R}\to\mathbb{R},t\mapsto f(x^{*}+tv)$
is differentiable at $t=0$ with derivative $0$. Then for
$t^{*}=\min\left\\{t>0\middle|\;H(x^{*}+tv)\neq H(x^{*})\right\\}$ it holds
that $f(x^{*})=f(x^{*}+tv)$ and
$\textnormal{df}^{C}(x^{*}+tv)<\textnormal{df}^{C}(x^{*})$. Repetition of this
procedure leads to a point $x^{**}$ with $\textnormal{df}^{C}(x^{**})=0$, i.e.
a vertex with $f(x^{**})=f(x^{*})$ after a finite number of steps. ∎

######  Lemma 28.

For $x\in\mathbb{R}^{n_{0}}$ and
$v\in\mathbb{R}^{n_{0}}\setminus\left\\{0\right\\}$ the limit $\lim_{t\searrow
0}\frac{1}{t}(f(x+tv)-f(x))$ exists. Furthermore assume that
$\varepsilon^{*}>0$ such that for all
$(l^{\prime}j^{\prime})\in\mathcal{I}\setminus C(x)$, $H(y)_{lj}$ is constant
for $y\in B_{\varepsilon^{*}}(x)$ and let $t^{*}>0$ small enough such that
$x+t^{*}v\in B_{\varepsilon^{*}}$. For $s\in\mathcal{S}^{C}(x+t^{*}v)$ it
holds that

$\lim_{t\searrow 0}\frac{1}{t}\left(f(x+tv)-f(x)\right)=\lim_{t\searrow
0}\frac{1}{t}\left(\tilde{f}_{s}(x+tv)-\tilde{f}_{s}(x)\right)=\langle
v,\nabla_{s}\rangle$

###### Proof.

Note that $H(x)_{lj}=H(x+t^{*}v)_{lj}$ for $(l,j)\in\mathcal{I}\setminus C(x)$
and $H(x)_{lj}=0$ for $(l,j)\in C(x)$. Hence $s$ is compatible with $H(x)$ and
by assumption $s$ is compatible with $H(x+t^{*}v)$. By Theorem 2 it follows
that $s$ is compatible with both $\tilde{H}_{s}(x)$ and
$\tilde{H}_{s}(x+t^{*}v)$. By equation (12) and the fact that all coordinates
of $y\mapsto\tilde{A}_{s}(y)$ are affine functions it follows that $s$ is
compatible with all $\tilde{H}_{s}(x+tv)$ for $t\in[0,t^{*}]$ and hence again
by Theorem 2 with all $H(x+tv)$ for $t\in[0,t^{*}]$. Applying the same theorem
again yields $A(x+tv)=\tilde{A}_{s}(x+tv)$ and hence
$f(x+tv)=\tilde{f}_{s}(x+tv)$ for all $t\in[0,t^{*}]$. The result now follows
from the fact that $\tilde{f}_{s}$ is affine with gradient exactly
$\nabla_{s}$ as in equation (34). ∎

######  Corollary 29.

Let $x\in\mathbb{R}^{n_{0}}$ be a regular vertex and $\varepsilon^{*}>0$ such
that for all $(l^{\prime}j^{\prime})\in\mathcal{I}\setminus C(x)$, $H(y)$ is
constant for $y\in B_{\varepsilon^{*}}(x)$. For every $v\in
B_{\varepsilon^{*}}(0)$ there exist non-negative coefficients
$\alpha_{1},\ldots,\alpha_{n_{0}}\in[0,\infty)$ and $n_{0}$ region separating
axes $a_{1},\ldots,a_{n_{0}}\in\mathcal{A}(x)$ such that

$f(x+v)=f(x)+\sum_{i=1}^{n_{0}}\alpha_{i}\lim_{t\searrow
0}\frac{1}{t}\left(f(x+ta_{i})-f(x)\right).$

###### Proof.

Let $v\in B_{\varepsilon^{*}}(0)$ and $s\in\mathcal{S}^{C}(x)$ such that $s$
is compatible with $H(x+v)$. Then $v\in\tilde{D}_{s}(x)$ by Corollary 27.
Lemma 17 and equation (31) imply that there exist coefficients
$\alpha_{lj}\in[0,\infty)$, $(l,j)\in C(x)$ such that

$v=\sum_{(l,j)\in\tilde{C}_{s}(x)}\alpha_{lj}\tilde{a}_{x,s,l,j}=\sum_{(l,j)\in\tilde{C}_{s}(x)}\alpha_{lj}a_{x,s_{lj},l,j}$
(55)

This implies that

$\displaystyle f(x+v)$
$\displaystyle\overset{(1)}{=}\tilde{f}_{s}(x+v)\overset{(2)}{=}\tilde{f}_{s}(x+\sum_{(l,j)\in\tilde{C}_{s}(x)}\alpha_{lj}a_{x,s_{lj},l,j})\overset{(3)}{=}\tilde{f}_{s}(x)+\langle\sum_{(l,j)\in\tilde{C}_{s}(x)}\alpha_{lj}a_{x,s_{lj},l,j},\nabla_{s}\rangle$
$\displaystyle\overset{}{=}f(x)+\sum_{(l,j)\in\tilde{C}_{s}(x)}\alpha_{lj}\langle
a_{x,s_{lj},l,j},\nabla_{s}\rangle\overset{(4)}{=}f(x)+\sum_{(l,j)\in\tilde{C}_{s}(x)}\alpha_{lj}\left(\tilde{f}_{s}(x+a_{x,s_{lj,l,j}})-\tilde{f}_{s}(x)\right)$
$\displaystyle\overset{(5)}{=}f(x)+\sum_{(l,j)\in\tilde{C}_{s}(x)}\alpha_{lj}\lim_{t\searrow
0}\frac{1}{t}\left(\tilde{f}_{s}(x+ta_{x,s_{lj},l,j})-\tilde{f}_{s}(x)\right)$
$\displaystyle\overset{(6)}{=}f(x)+\sum_{(l,j)\in\tilde{C}_{s}(x)}\alpha_{lj}\lim_{t\searrow
0}\frac{1}{t}\left(f(x+ta_{x,s_{lj},l,j})-f(x)\right).$

The justification of the steps is as follows: (1) follows from Corollary 3,
(2) from equation (55), (3), (4) and (5) from the fact that the functions
$f_{s^{\prime}}$, $s^{\prime}\in\mathcal{S}$ are affine and (6) is an
implication of Lemma 28. Furthermore $|\tilde{C}_{s}(x)|=n_{0}$ by the
regularity assumption. ∎

###### Proof of Proposition 20.

The claim follows from Lemmas 28 and Corollary 29. ∎

##  Appendix B Bounds on the number of affine regions

In the sequel, let $N(n_{0},\dots,n_{L})$ denote the maximal number of affine
regions that a ReLU feed-forward neural network of depth $L$ and widths
$n_{0},\dots,n_{L}$ can have. We will now present some bounds on
$N(n_{0},\dots,n_{L})$.

The work [1] establishes the following bound on the number of affine regions:

######  Lemma 30 (Montúfar).

Let $L\in\mathbb{N}$. A neural network with $L+2$ layers ($L$ hidden layers)
of widths $n_{0},\dots,n_{L+1}$ is a piece-wise affine function with at most

$N(n_{0},\dots,n_{L})\leq\prod_{i=1}^{L}\sum_{j=0}^{\min(n_{0},\dots,n_{i})}{n_{i}\choose
j}$ (56)

affine regions.

It uses a result from 1943 by R.C.Buck, which bounds the number of regions
that are carved out of $\mathbb{R}^{n_{0}}$ by $n_{1}$ hyperplanes, see [13]:

######  Lemma 31 (Buck).

In $\mathbb{R}^{d}$, $N\in\mathbb{N}$ hyperplanes cut the space
$\mathbb{R}^{d}$ into at most $\sum_{j=0}^{d}{N\choose j}$ regions.

This bound was improved in [2] to

$N(n_{0},\dots,n_{L})\leq\sum_{j_{1},\dots,j_{L}\in
J}\prod_{i=1}^{L}{n_{i}\choose j},$ (57)

where $J=\left\\{(j_{1},\dots,j_{L})\in\mathbb{N}^{L}|\;\forall
i\in\left\\{1,\dots,L\right\\}:\;j_{i}\leq\min(n_{0},n_{1}-j_{1},\dots,n_{L-1}-j_{L-1},n_{L})\right\\}$.
Later, an abstract framework for the construction of such upper bounds was
presented in [3]. It yields bounds formulated as the L1 norm of a product of
matrices

$N(n_{0},\dots,n_{L})\leq\|B^{(\gamma)}_{n_{L}}M_{n_{L-1},n_{L}}\dots
B^{(\gamma)}_{n_{1}}M_{n_{0},n_{1}}e^{(n_{0}+1)}_{n_{0}+1}\|_{1},$ (58)

where $M_{n_{0},n_{1}},\dots,M_{n_{L-1},n_{L}}$ are defined by
$M_{a,b}=\left(\delta_{i,\min(j,a+1)}\right)_{i\in\left\\{1,\dots,b+1\right\\},j\in\left\\{1,\dots,a+1\right\\}}\in\mathbb{N}^{b+1\times
a+1}$ for $a,b\in\mathbb{N}_{+}$,
$e^{(n_{0}+1)}_{n_{0}+1}=\left(\delta_{n_{0}+1,i}\right)_{i\in\\{1,\dots,n_{0}+1\\}}\in\mathbb{N}^{n_{0}+1}$
and $B^{(\gamma)}_{a}\in\mathbb{N}^{a+1\times a+1}$ are upper triangular
matrices depending on some family of vectors
$(\gamma_{ij})_{i,j\in\mathbb{N}}$. The bounds (56) and (57) can be derived
from (58), see [3].

## References

  * [1] G. Montúfar, R. Pascanu, K. Cho, and Y. Bengio, “On the number of linear regions of deep neural networks,” in Proceedings of the 27th International Conference on Neural Information Processing Systems - Volume 2, NIPS’14, (Cambridge, MA, USA), pp. 2924–2932, MIT Press, 2014.
  * [2] T. Serra, C. Tjandraatmadja, and S. Ramalingam, “Bounding and counting linear regions of deep neural networks,” CoRR, vol. abs/1711.02114, 2017.
  * [3] P. Hinz and S. van de Geer, “A framework for the construction of upper bounds on the number of affine linear regions of relu feed-forward neural networks,” arXiv preprint arXiv:1806.01918, 2018.
  * [4] R. Fourer, “A simplex algorithm for piecewise-linear programming i: Derivation and proof,” Mathematical Programming, vol. 33, pp. 204–233, Nov 1985.
  * [5] R. Koenker and G. Bassett Jr, “Regression quantiles,” Econometrica: journal of the Econometric Society, pp. 33–50, 1978.
  * [6] P. Bloomfield and W. L. Steiger, Least Absolute Deviations: Theory, Applications and Algorithms, pp. 37–76.  Boston, MA: Birkhäuser Boston, 1983.
  * [7] I. Barrodale and F. D. Roberts, “An improved algorithm for discrete l1 linear approximation,” SIAM Journal on Numerical Analysis, vol. 10, no. 5, pp. 839–848, 1973.
  * [8] J. L. Powell, “Least absolute deviations estimation for the censored regression model,” Journal of Econometrics, vol. 25, no. 3, pp. 303–325, 1984.
  * [9] M. Buchinsky, “Changes in the u.s. wage structure 1963-1987: Application of quantile regression,” Econometrica, vol. 62, no. 2, pp. 405–458, 1994\.
  * [10] B. Fitzenberger, “Computational aspects of censored quantile regression,” Lecture Notes-Monograph Series, vol. 31, pp. 171–186, 1997.
  * [11] R. Tibshirani, “Regression shrinkage and selection via the lasso,” Journal of the Royal Statistical Society: Series B (Methodological), vol. 58, no. 1, pp. 267–288, 1996.
  * [12] N. Qian, “On the momentum term in gradient descent learning algorithms,” Neural networks, vol. 12, no. 1, pp. 145–151, 1999.
  * [13] R. C. Buck, “Partition of space,” American Mathematical Monthly, vol. 50, no. 9, pp. 541–544, 1943.


# Complex Dynamics of $\displaystyle{z_{n+1}=\frac{\alpha+\beta z_{n}+\gamma
z_{n-1}}{A+Bz_{n}+Cz_{n-1}}}$

Sk. Sarif Hassan
International Centre for Theoretical Sciences
Tata Institute of Fundamental Research
Bangalore $560012$, India
Email<EMAIL_ADDRESS>

###### Abstract

The dynamics of the second order rational difference equation
$\displaystyle{z_{n+1}=\frac{\alpha+\beta z_{n}+\gamma
z_{n-1}}{A+Bz_{n}+Cz_{n-1}}}$ with complex parameters and arbitrary complex
initial conditions is investigated. In the complex set up, the local
asymptotic stability and boundedness are studied vividly for this difference
equation. Several interesting characteristics of the solutions of this
equation, using computations, which does not arise when we consider the same
equation with positive real parameters and initial conditions are shown. The
chaotic solutions of the difference equation is absolutely new feature in the
complex set up which is also shown in this article. Some of the interesting
observations led us to pose some open interesting problems regarding chaotic
and higher order periodic solutions and global asymptotic convergence of this
equation.

Keywords: Rational difference equation, Local asymptotic stability, Chaotic
trajectory and Periodicity.
Mathematics Subject Classification: 39A10, 39A11

##  1 Introduction and Preliminaries

A rational difference equation is a nonlinear difference equation of the form

$z_{n+1}=\frac{\alpha+\sum_{i=0}^{k}\beta_{i}z_{n-i}}{A+\sum_{i=0}^{k}B_{i}z_{n-i}}$

where the initial conditions $z_{0},z_{-1},\dots,z_{-k}$ are such that the
denominator never vanishes for any $n$.
Consider the equation

$z_{n+1}=\frac{\alpha+\beta z_{n}+\gamma
z_{n-1}}{A+Bz_{n}+Cz_{n-1}},n=0,1,\ldots$ (1)

where all the parameters and the initial conditions $z_{-1}$ and $z_{0}$ are
arbitrary complex number.

This second order rational difference equation Eq.(1) is studied when the
parameters are real numbers and initial conditions are non-negative real
numbers in [2]. In this present article it is an attempt to understand the
same in the complex plane.

Here, a very brief review of some well known results which will be useful in
order to apprehend the behavior of solutions of the difference equation (1).

Let $f:\mathbb{D}^{2}\rightarrow\mathbb{D}$ where
$\mathbb{D}\subseteq\mathbf{C}$ be a continuously differentiable function.
Then for any pair of initial conditions $z_{0},z_{-1}\in\mathbb{D}$, the
difference equation

$z_{n+1}=f(z_{n},z_{n-1})\hskip 18.06749pt\hskip 18.06749pt$ (2)

with initial conditions $z_{-1},z_{0}\in\mathbb{D}.$

Then for any _initial value_ , the difference equation (1) will have a unique
solution $\\{z_{n}\\}_{n}$.

A point $\overline{z}$ $\in$ $\mathbb{D}$ is called equilibrium point of
Eq.(2) if

$f(\overline{z},\overline{z})=\overline{z}.$

The _linearized equation_ of Eq.(2) about the equilibrium $\bar{z}$ is the
linear difference equation

$\displaystyle{z_{n+1}=a_{0}z_{n}+a_{1}z_{n-1}\hskip 18.06749pt,\hskip
18.06749ptn=0,1,\ldots}$ (3)

where for $i=0$ and $1$.

$a_{i}=\frac{\partial f}{\partial u_{i}}(\overline{z},\overline{z}).$

The _characteristic equation_ of Eq.(2) is the equation

$\lambda^{2}-a_{0}\lambda-a_{1}=0.$ (4)

The following are the briefings of the linearized stability criterions which
are useful in determining the local stability character of the equilibrium
$\overline{z}$ of Eq.(2), [1].

Let $\overline{z}$ be an equilibrium of the difference equation
$z_{n+1}=f(z_{n},z_{n-1})$.

  * •

The equilibrium $\bar{z}$ of Eq. (2) is called locally stable if for every
$\epsilon>0$, there exists a $\delta>0$ such that for every $z_{0}$ and
$z_{-1}$ $\in$ $\mathbb{C}$ with $|z_{0}-\bar{z}|+|z_{-1}-\bar{z}|<\delta$ we
have $|z_{n}-\bar{z}|<\epsilon$ for all $n>-1$.

  * •

The equilibrium $\bar{z}$ of Eq. (2) is called locally stable if it is locally
stable and if there exist a $\gamma>0$ such that for every $z_{0}$ and
$z_{-1}$ $\in$ $\mathbb{C}$ with $|z_{0}-\bar{z}|+|z_{-1}-\bar{z}|<\gamma$ we
have $\lim_{n\to\infty}z_{n}=\bar{z}$.

  * •

The equilibrium $\bar{z}$ of Eq. (2) is called global attractor if for every
$z_{0}$ and $z_{-1}$ $\in$ $\mathbb{C}$, we have
$\lim_{n\to\infty}z_{n}=\bar{z}$.

  * •

The equilibrium of equation Eq. (2) is called globally asymptotically
stable/fit is stable and is a global attractor.

  * •

The equilibrium $\bar{z}$ of Eq. (2) is called unstable if it is not stable.

  * •

The equilibrium $\bar{z}$ of Eq. (2) is called source or repeller if there
exists $r>0$ such that for every $z_{0}$ and $z_{-1}$ $\in$ $\mathbb{C}$ with
$|z_{0}-\bar{z}|+|z_{-1}-\bar{z}|<r$ we have $|z_{n}-\bar{z}|\geq r$. Clearly
a source is an unstable equilibrium.

Result 1.1: (Clark’s Theorem) The sufficient condition for the asymptotic
stability of the difference equation (1) is

$\left\lvert a_{0}\right\rvert+\left\lvert a_{1}\right\rvert<1$

##  2 Difference Equation and Its Transformed Forms

The following difference equation is considered to be studied here.

$z_{n+1}=\frac{\alpha+\beta z_{n}+\gamma
z_{n-1}}{A+Bz_{n}+Cz_{n-1}},n=0,1,\ldots$ (5)

where all the parameters are complex number and the initial conditions
$z_{-1}$ and $z_{0}$ are arbitrary complex numbers.

We will consider three different cases of the Eq.(2) which are as follows:

###  2.1 The case $\beta=\gamma=0:z_{n+1}=\frac{\alpha}{A+Bz_{n}+Cz_{n-1}}$

By the change of variables, $z_{n}=\frac{\alpha}{A}w_{n}$, the difference
equation $z_{n+1}=\frac{\alpha}{A+Bz_{n}+Cz_{n-1}}$ reduced to the difference
equation

$w_{n+1}=\frac{1}{1+pw_{n}+qw_{n-1}},n=0,1,\ldots$ (6)

where $p=\frac{\alpha B}{A^{2}}$ and $q=\frac{\alpha C}{A^{2}}$.

###  2.2 The case $\alpha=\gamma=0:z_{n+1}=\frac{\beta
z_{n}}{A+Bz_{n}+Cz_{n-1}}$

By the change of variables, $z_{n}=\frac{A}{C}w_{n}$, the difference equation
$z_{n+1}=\frac{\beta z_{n}}{A+Bz_{n}+Cz_{n-1}}$ reduced to the difference
equation

$w_{n+1}=\frac{w_{n}}{1+pw_{n}+qw_{n-1}},n=0,1,\ldots$ (7)

where $p=\frac{\beta}{A}$ and $q=\frac{B}{C}$.

###  2.3 The case $\alpha=\beta=0:z_{n+1}=\frac{\gamma
z_{n-1}}{A+Bz_{n}+Cz_{n-1}}$

By the change of variables, $z_{n}=\frac{\gamma}{C}w_{n}$, the difference
equation $z_{n+1}=\frac{\gamma z_{n-1}}{A+Bz_{n}+Cz_{n-1}}$ reduced to the
difference equation

$w_{n+1}=\frac{w_{n-1}}{p+qw_{n}+w_{n-1}},n=0,1,\ldots$ (8)

where $p=\frac{A}{\gamma}$ and $q=\frac{B}{C}$.

Without any loss of generality, we shall now onward focus only on the three
difference equations (6), (7) and (8).

##  3 Local Asymptotic Stability of the Equilibriums

In this section we establish the local stability character of the equilibria
of Eq.(1) in three difference cases as stated in the section 2.

###  3.1 Local Asymptotic Stability of $w_{n+1}=\frac{1}{1+pw_{n}+qw_{n-1}}$

The equilibrium points of Eq.(6) are the solutions of the quadratic equation

$\bar{w}=\frac{1}{1+p\bar{w}+q\bar{w}}$

Eq.(6) has the two equilibria points
$\bar{w}_{1,2}=\frac{-1-\sqrt{1+4p+4q}}{2(p+q)}$ and
$\frac{-1+\sqrt{1+4p+4q}}{2(p+q)}$ respectively. The linearized equation of
the rational difference equation(6) with respect to the equilibrium point
$\bar{w}_{1}=\frac{-1-\sqrt{1+4p+4q}}{2(p+q)}$ is

$\displaystyle{z_{n+1}+\frac{4p}{\left(-1+\sqrt{1+4p+4q}\right)^{2}}w_{n}+\frac{4q}{\left(-1+\sqrt{1+4p+4q}\right)^{2}}w_{n-1}=0,n=0,1,\ldots}$
(9)

with associated characteristic equation

$\lambda^{2}+\frac{4p}{\left(-1+\sqrt{1+4p+4q}\right)^{2}}\lambda+\frac{4q}{\left(-1+\sqrt{1+4p+4q}\right)^{2}}=0.$
(10)

The following result gives the local asymptotic stability of the equilibrium
$\bar{w}_{1}$ of the Eq. (6).

######  Theorem 3.1.

The equilibriums $\bar{w}_{1}=\frac{-1-\sqrt{1+4p+4q}}{2(p+q)}$ of Eq.(6) is

locally asymptotically stable if

$\left\lvert\frac{4p}{\left(-1+\sqrt{1+4p+4q}\right)^{2}}\right\rvert+\left\lvert{\frac{4q}{\left(-1+\sqrt{1+4p+4q}\right)^{2}}}\right\rvert<1$

###### Proof.

The zeros of the characteristic equation (10) has two zeros which are
$-\frac{2q}{p+\sqrt{p^{2}-4pq+2q\left(-1-2q+\sqrt{1+4p+4q}\right)}}$ and
$\frac{2q}{-p+\sqrt{p^{2}-4pq+2q\left(-1-2q+\sqrt{1+4p+4q}\right)}}$.
Therefore by _Clark’s theorem_ , the equilibrium
$\bar{w}_{1}=\frac{-1-\sqrt{1+4p+4q}}{2(p+q)}$ is _locally asymptotically
stable_ if the sum of the modulus of two coefficients is less than $1$.
Therefore the condition of the polynomial (10) reduces to
$\left\lvert\frac{4p}{\left(-1+\sqrt{1+4p+4q}\right)^{2}}\right\rvert+\left\lvert{\frac{4q}{\left(-1+\sqrt{1+4p+4q}\right)^{2}}}\right\rvert<1$.

∎

The linearized equation of the rational difference equation (6) with respect
to the equilibrium point $\bar{w}_{2}=\frac{-1+\sqrt{1+4p+4q}}{2(p+q)}$ is

$\displaystyle{w_{n+1}+\frac{4p}{\left(1+\sqrt{1+4p+4q}\right)^{2}}w_{n}+\frac{4q}{\left(1+\sqrt{1+4p+4q}\right)^{2}}w_{n-1}=0,n=0,1,\ldots}$
(11)

with associated characteristic equation

$\lambda^{2}+\frac{4p}{\left(1+\sqrt{1+4p+4q}\right)^{2}}\lambda+-\frac{4q}{\left(1+\sqrt{1+4p+4q}\right)^{2}}=0.$
(12)

######  Theorem 3.2.

The equilibriums $\bar{w}_{2}=\frac{-1+\sqrt{1+4p+4q}}{2(p+q)}$ of Eq.(6) is

locally asymptotically stable if

$\left\lvert\frac{4p}{\left(1+\sqrt{1+4p+4q}\right)^{2}}\right\rvert+\left\lvert\frac{4q}{\left(1+\sqrt{1+4p+4q}\right)^{2}}\right\rvert<1$

###### Proof.

Proof the theorem follows from _Clark’s theorem_ of local asymptotic stability
of the equilibriums. The condition for the local asymptotic stability reduces
to
$\left\lvert\frac{4p}{\left(1+\sqrt{1+4p+4q}\right)^{2}}\right\rvert+\left\lvert\frac{4q}{\left(1+\sqrt{1+4p+4q}\right)^{2}}\right\rvert<1$.

∎

Here is an example case for the local asymptotic stability of the
equilibriums.

For $p=\frac{1}{2}$ and $q=\frac{i}{2}$ the equilibriums are
$-1.6838+1.13355i$ and $0.683802-0.133552i$. For the equilibrium
$1.6838+1.13355i$, the coefficients of the characteristic polynomial (10) are
$0.775125+1.90868i$ and $-1.90868-0.775125i$ with same modulus $2.06006$.
Therefore the condition as stated in the _Theorem 3.1_ does not hold.
Therefore the equilibrium $1.6838+1.13355i$ is _unstable_.

For the equilibrium $0.683802-0.133552i$, the coefficients of the
characteristic polynomial (12) are $-0.224875+0.091323i$ and
$-0.091323-0.224875i$ with same modulus $0.242711$. Therefore the condition as
stated in the _Theorem 3.2_ is hold good. Therefore the equilibrium
$0.683802-0.133552i$ is _locally asymptotically stable_.

It is seen that in the case of real positive parameters and initials values,
the positive equilibrium of the difference equation (6) is globally
asymptotically stable [2]. But the result is not holding well in the complex
set up.

###  3.2 Local Asymptotic Stability of
$w_{n+1}=\frac{w_{n}}{1+pw_{n}+qw_{n-1}}$

The equilibrium points of Eq.(7) are the solutions of the quadratic equation

$\bar{w}=\frac{\bar{w}}{1+p\bar{w}+q\bar{w}}$

The Eq.(7) has only the zero equilibrium. The linearized equation of the
rational difference equation(7) with respect to the zero equilibrium is

$\displaystyle{w_{n+1}=w_{n},n=0,1,\ldots}$ (13)

with associated characteristic equation

$\lambda^{2}-\lambda=0.$ (14)

The following result gives the local asymptotic stability of the zero
equilibrium of the Eq. (7).

######  Theorem 3.3.

The zero equilibriums of the Eq. (7) is non-hyperbolic.

###### Proof.

The zeros of the characteristic equation (14) has two zeros which are $0$ and
$1$. Therefore by definition, the zero equilibrium is non-hyperbolic as the
modulus of one zero is $1$. ∎

It is nice to note that in the case of real positive parameters and initials
values, the zero equilibrium of the difference equation (7) is globally
asymptotically stable for the parameter $p\geq 1$ [2]. But in the case of
complex, the zero equilibrium is non-hyperbolic as we have seen the previous
theorem.

###  3.3 Local Asymptotic Stability of
$w_{n+1}=\frac{w_{n-1}}{p+qw_{n}+w_{n-1}}$

The equilibrium points of Eq.(8) are the solutions of the quadratic equation

$\bar{w}=\frac{\bar{w}}{p+q\bar{w}+\bar{w}}$

The Eq.(8) has two equilibriums which are $0$ and $\frac{1-p}{1+q}$. The
linearized equation of the rational difference equation(8) with respect to the
zero equilibrium is

$\displaystyle{w_{n+1}=\frac{1}{p}w_{n},n=0,1,\ldots}$ (15)

with associated characteristic equation

$\lambda^{2}-\frac{\lambda}{p}=0.$ (16)

The following result gives the local asymptotic stability of the zero
equilibrium of the Eq. (8).

######  Theorem 3.4.

The zero equilibriums of the Eq. (8) is locally asymptotically stable if
$\left\lvert p\right\rvert\geq 1$ and repeller if $\left\lvert
p\right\rvert<1$.

###### Proof.

The zeros of the characteristic equation (15) has two zeros which are $0$ and
$\frac{1}{p}$. Therefore by definition, the zero equilibrium is locally
asymptotically stable if $\left\lvert\frac{1}{p}\right\rvert<1$ and unstable
(repeller) if $\left\lvert\frac{1}{p}\right\rvert\geq 1$. Hence the required
is followed. ∎

The linearized equation of the rational difference equation(8) with respect to
the equilibrium $\frac{1-p}{1+q}$ is

$\displaystyle{w_{n+1}=\frac{1+pq}{1+q}w_{n}+\frac{p-1}{1+q}w_{n-1}=0,n=0,1,\ldots}$
(17)

with associated characteristic equation

$\lambda^{2}-\frac{1+pq}{1+q}\lambda-\frac{p-1}{1+q}=0.$ (18)

The following result gives the local asymptotic stability of the equilibrium
$\frac{1-p}{1+q}$ of the Eq. (8).

######  Theorem 3.5.

The zero equilibriums of the Eq. (8) is locally asymptotically stable if

$\left\lvert 1+pq\right\rvert+\left\lvert p-1\right\rvert<\left\lvert
1+q\right\rvert$

###### Proof.

The equilibrium $\frac{1-p}{1+q}$ of the characteristic equation (18) would be
_locally asymptotically stable_ if the sum of the modulus of the coefficients
of the characteristic equation (18) is less than 1. That is by _Clark’s
theorem_ ,
$\left\lvert\frac{1+pq}{1+q}\right\rvert+\left\lvert\frac{p-1}{1+q}\right\rvert<1$,
that is

$\left\lvert 1+pq\right\rvert+\left\lvert p-1\right\rvert<\left\lvert
1+q\right\rvert$

∎

Here is an example case for the local asymptotic stability of the
equilibriums.

For $p=1+\frac{i}{2}$ $(\left\lvert p\right\rvert=1.1180>1$) and
$q=\frac{1}{10}+i$ ($\left\lvert q\right\rvert=1.005>1$) the equilibriums are
$0$ and $-0.226244-0.248869i$. For the equilibrium $-0.226244-0.248869i$, the
coefficients of the characteristic polynomial (18) are $0.533622+0.545503i$
and $-0.203409-0.0330197i$ with modulus $0.763103$ and $0.206072$
respectively. Therefore the condition as stated in the _Theorem 3.5_ hold
good. Therefore the equilibrium $1.6838+1.13355i$ is _locally asymptotically
stable_.

In the case of real positive parameters and intimal values of the difference
equation (8), the positive equilibrium is locally asymptotically stable if
$p<1$ and $q<1$ but in the complex set, it is encountered through the example
above is that the equilibrium $\frac{1-p}{1+q}$ is locally asymptotically
stable even though $\left\lvert p\right\rvert>1$ and $\left\lvert
q\right\rvert>1$ [2].

##  4 Boundedness

In this section we would like to explore the boundedness of the solutions of
the three difference equations (6), (7) and (8).

Now we would like to try to find open ball $B(0,\epsilon)\in\mathbf{C}$ such
that if $w_{n}\in B(0,\epsilon)$ and $w_{n-1}\in B(0,\epsilon)$ then
$w_{n+1}\in B(0,\epsilon)$ for all $n\geq 0$. In other words, if the initial
values $w_{0}$ and $w_{-1}$ belong to $B(0,\epsilon)$ then the solution
generated by the difference equations would essentially be within the open
ball $B(0,\epsilon)$.

######  Theorem 4.1.

For the difference equation (6), for every $\epsilon>0$, if $w_{n}$ and
$w_{n-1}$ $\in B(0,\epsilon)$ then $w_{n+1}\in B(0,\epsilon)$ provided

$\left\lvert p\right\rvert\geq 1+\left\lvert q\right\rvert$

###### Proof.

Let $\\{w_{n}\\}$ be a solution of the equation Eq.(6). Let $\epsilon>0$ be
any arbitrary real number. Consider $w_{n},w_{n-1}\in B(0,\epsilon)$. We need
to find out an $\epsilon$ such that $w_{n+1}\in B(0,\epsilon)$ for all $n$. It
is follows from the Eq.(6) that for any $\epsilon>0$, using Triangle
inequality for

$\left\lvert
w_{n+1}\right\rvert=\left\lvert\frac{1}{1+pw_{n}+qw_{n-1}}\right\rvert\leq\left\lvert\frac{1}{pw_{n}+qw_{n-1}}\right\rvert\leq\frac{1}{(\left\lvert
p\right\rvert-\left\lvert q\right\rvert)\epsilon}$

In order to ensure that $\left\lvert w_{n+1}\right\rvert<\epsilon^{\prime}$,
(Assuming $\epsilon^{\prime}=\frac{1}{\epsilon}$) it is needed to be

$\frac{1}{(\left\lvert p\right\rvert-\left\lvert q\right\rvert)}<1$

That is $\left\lvert p\right\rvert\geq 1+\left\lvert q\right\rvert$. Therefore
the required is followed.

∎

######  Theorem 4.2.

For the difference equation (7), for every $\epsilon>0$, if $w_{n}$ and
$w_{n-1}$ $\in B(0,\epsilon)$ then $w_{n+1}\in B(0,\epsilon)$ provided

$\left\lvert p\right\rvert\geq\left\lvert q\right\rvert+\frac{1}{\epsilon}$

###### Proof.

Let $\\{w_{n}\\}$ be a solution of the equation Eq.(7). Let $\epsilon>0$ be
any arbitrary real number. Consider $w_{n},w_{n-1}\in B(0,\epsilon)$. We need
to find out an $\epsilon$ such that $w_{n+1}\in B(0,\epsilon)$ for all $n$. It
is follows from the Eq.(7), that for any $\epsilon>0$, using Triangle
inequality for

$\left\lvert
w_{n+1}\right\rvert=\left\lvert\frac{w_{n}}{1+pw_{n}+qw_{n-1}}\right\rvert\leq\frac{\epsilon}{(\left\lvert
p\right\rvert-\left\lvert q\right\rvert)\epsilon}$

Therefore, $\left\lvert w_{n+1}\right\rvert<\frac{1}{\left\lvert
p\right\rvert-\left\lvert q\right\rvert}$. We need $\frac{1}{\left\lvert
p\right\rvert-\left\lvert q\right\rvert}$ to be less than $\epsilon$.
Therefore,

$\left\lvert p\right\rvert\geq\left\lvert q\right\rvert+\frac{1}{\epsilon}$

is followed.

∎

######  Theorem 4.3.

For the difference equation (8), for every $\epsilon>0$, if $w_{n}$ and
$w_{n-1}$ $\in B(0,\epsilon)$ then $w_{n+1}\in B(0,\epsilon)$ and $\left\lvert
p\right\rvert<1$ provided

$\epsilon<\frac{\left\lvert p\right\rvert-1}{\left\lvert q\right\rvert+1}$

###### Proof.

Let $\\{w_{n}\\}$ be a solution of the equation Eq.(8). Let $\epsilon>0$ be
any arbitrary real number and $\left\lvert p\right\rvert>1$. Consider
$w_{n},w_{n-1}\in B(0,\epsilon)$. We need to find out an $\epsilon$ such that
$w_{n+1}\in B(0,\epsilon)$ for all $n$. It is follows from the Eq.(8), that
for any $\epsilon>0$, using Triangle inequality for

$\left\lvert
w_{n+1}\right\rvert=\left\lvert\frac{w_{n-1}}{p+qw_{n}+w_{n-1}}\right\rvert\leq\frac{\epsilon}{\left\lvert
p\right\rvert-\left\lvert q\right\rvert\epsilon-\epsilon}$

Therefore,

$\left\lvert w_{n+1}\right\rvert\leq\frac{\epsilon}{\left\lvert
p\right\rvert-\left\lvert q\right\rvert\epsilon-\epsilon}$

In order to ensure that $\left\lvert w_{n+1}\right\rvert<\epsilon$, it is
needed to be

$\frac{1}{\left\lvert p\right\rvert-\left\lvert
q\right\rvert\epsilon-\epsilon}<1$

That is

$\epsilon<\frac{\left\lvert p\right\rvert-1}{\left\lvert q\right\rvert+1}$

Therefore the required is followed.

∎

##  5 Periodic of Solutions

A solution $\\{w_{n}\\}_{n}$ of a difference equation is said to be _globally
periodic_ of period $t$ if $w_{n+t}=w_{n}$ for any given initial conditions.
solution $\\{w_{n}\\}_{n}$ is said to be _periodic with prime period_ $p$ if p
is the smallest positive integer having this property.

We shall first look for the prime period two solutions of the three difference
equations (6), (7) and (8) and their corresponding local stability analysis.

###  5.1 Prime Period Two Solutions of Eq. (6)

Let $\ldots,\phi,~{}\psi,~{}\phi,~{}\psi,\ldots$, $\phi\neq\psi$ be a prime
period two solution of the difference equation
$w_{n+1}=\frac{1}{1+pw_{n}+qw_{n-1}}$. Then $\phi=\frac{1}{1+p\psi+q\phi}$ and
$\psi=\frac{1}{1+p\phi+q\psi}$. This two equations lead to the set of
solutions (prime period two) except the equilibriums as
$\left\\{\phi\to\frac{0.5\left(p-q+\sqrt{p^{2}+p(-2-4q)q+q^{2}+4q^{3}}\right)}{q(-p+q)},\psi\to\frac{0.5p-0.5q-0.5\sqrt{p^{2}+p(-2-4q)q+q^{2}+4q^{3}}}{q(-p+q)}\right\\}$
and
$\left\\{\phi\to\frac{q}{-0.5p+0.5q-0.5\sqrt{p^{2}+p(-2-4q)q+q^{2}+4q^{3}}},\psi\to\frac{-0.5-\frac{0.5\sqrt{p^{2}+p(-2-4q)q+q^{2}+4q^{3}}}{p-q}}{q}\right\\}$.

Let $\ldots,\phi,~{}\psi,~{}\phi,~{}\psi,\ldots$, $\phi\neq\psi$ be a prime
period two solution of the equation (6). We set

$u_{n}=w_{n-1}$ $v_{n}=w_{n}$

Then the equivalent form of the difference equation (6) is

$u_{n+1}=v_{n}$ $v_{n+1}=\frac{1}{1+pv_{n}+qu_{n}}$

Let T be the map on $\mathbb{C}\times\mathbb{C}$ to itself defined by

$T\left(\begin{array}[]{c}u\\\ v\\\
\end{array}\right)=\left(\begin{array}[]{c}v\\\ \frac{1}{1+pv+qu}\\\
\end{array}\right)$

Then $\left(\begin{array}[]{c}\phi\\\ \psi\\\ \end{array}\right)$ is a fixed
point of $T^{2}$, the second iterate of $T$.

$T^{2}\left(\begin{array}[]{c}u\\\ v\\\
\end{array}\right)=\left(\begin{array}[]{c}\frac{1}{1+pv+qu}\\\ \\\
\frac{1}{1+p\frac{1}{1+pv+qu}+qv}\\\ \end{array}\right)$

$T^{2}\left(\begin{array}[]{c}u\\\ v\\\
\end{array}\right)==\left(\begin{array}[]{c}g(u,v)\\\ h(u,v)\\\
\end{array}\right)$

where $g(u,v)=\frac{1}{1+pv+qu}$ and
$h(u,v)=\frac{1}{1+\frac{p}{1+\text{pv}+\text{qu}}+\text{qv}}$. Clearly the
two cycle is locally asymptotically stable when the eigenvalues of the
Jacobian matrix $J_{T^{2}}$, evaluated at $\left(\begin{array}[]{c}\phi\\\
\psi\\\ \end{array}\right)$ lie inside the unit disk.

We have,

$J_{T^{2}}\left(\begin{array}[]{c}\phi\\\ \psi\\\
\end{array}\right)=\left(\begin{array}[]{cc}\frac{\delta g}{\delta
u}(\phi,\psi)&\frac{\delta g}{\delta v}(\phi,\psi)\\\ \\\ \frac{\delta
h}{\delta u}(\phi,\psi)&\frac{\delta h}{\delta v}(\phi,\psi)\\\
\end{array}\right)$

where $\frac{\delta g}{\delta u}(\phi,\psi)=-\frac{q}{(1+q\phi+p\psi)^{2}}$
and $\frac{\delta g}{\delta v}(\phi,\psi)=-\frac{p}{(1+q\phi+p\psi)^{2}}$

$\frac{\delta h}{\delta
u}(\phi,\psi)=\frac{pq}{\left((1+q\phi)(1+q\psi)+p\left(1+\psi+q\psi^{2}\right)\right)^{2}}$
and $\frac{\delta h}{\delta
v}(\phi,\psi)=\frac{-q+\frac{p^{2}}{(1+q\phi+p\psi)^{2}}}{\left(1+q\psi+\frac{p}{1+q\phi+p\psi}\right)^{2}}$

Now, set

$\chi=\frac{\delta g}{\delta u}(\phi,\psi)+\frac{\delta h}{\delta
v}(\phi,\psi)=-\frac{q}{(1+q\phi+p\psi)^{2}}+\frac{-q+\frac{p^{2}}{(1+q\phi+p\psi)^{2}}}{\left(1+q\psi+\frac{p}{1+q\phi+p\psi}\right)^{2}}$
$\lambda=\frac{\delta g}{\delta u}(\phi,\psi)\frac{\delta h}{\delta
v}(\phi,\psi)-\frac{\delta g}{\delta v}(\phi,\psi)\frac{\delta h}{\delta
u}(\phi,\psi)=\frac{q^{2}}{\left((1+q\phi)(1+q\psi)+p\left(1+\psi+q\psi^{2}\right)\right)^{2}}$

In particular for the prime period $2$ solution,

$\left\\{\phi\to\frac{0.5\left(p-q+\sqrt{p^{2}+p(-2-4q)q+q^{2}+4q^{3}}\right)}{q(-p+q)},\psi\to\frac{0.5p-0.5q-0.5\sqrt{p^{2}+p(-2-4q)q+q^{2}+4q^{3}}}{q(-p+q)}\right\\}$,
we shall see the local asymptotic stability for some example cases of
parameters $p$ and $q$. The general form of $\chi$ and $\lambda$ would be very
complected. Consider the prime period two solution of the difference equation
(6), $\phi\to-0.0843748+0.0622145i,\psi\to-0.0822456-0.0594374i$ corresponding
two the parameters $p\to 100+i$ and $q\to 6+0.1i$.

In this case, $\left\lvert\chi\right\rvert=0.0735211$ and
$\left\lvert\lambda\right\rvert=0.004075$. Therefore, by the Linear Stability
theorem ($|\chi|<1+|\lambda|<2$) the prime period $2$ solution
$\phi\to-0.0843748+0.0622145i,\psi\to-0.0822456-0.0594374i$ is _locally
asymptotically stable_.

###  5.2 Prime Period Two Solutions of Eq. (7)

Let $\ldots,\phi,~{}\psi,~{}\phi,~{}\psi,\ldots$, $\phi\neq\psi$ be a prime
period two solution of the difference equation
$w_{n+1}=\frac{w_{n}}{1+pw_{n}+qw_{n-1}}$. Then
$\phi=\frac{\psi}{1+p\psi+q\phi}$ and $\psi=\frac{\phi}{1+p\phi+q\psi}$. This
two equations lead to the set of solutions (prime period two) except the
equilibriums as
$\left\\{\phi\to\frac{1}{0.5p-0.5q-0.5\sqrt{p^{2}-q^{2}}},\psi\to-\frac{1}{q}+\frac{\sqrt{p^{2}-q^{2}}}{pq-q^{2}}\right\\},\left\\{\phi\to\frac{1}{0.5p-0.5q+0.5\sqrt{p^{2}-q^{2}}},\psi\to-\frac{1}{q}-\frac{\sqrt{p^{2}-q^{2}}}{pq-q^{2}}\right\\}$.

Let $\ldots,\phi,~{}\psi,~{}\phi,~{}\psi,\ldots$, $\phi\neq\psi$ be a prime
period two solution of the equation (7). We set

$u_{n}=w_{n-1}$ $v_{n}=w_{n}$

Then the equivalent form of the difference equation (7) is

$u_{n+1}=v_{n}$ $v_{n+1}=\frac{v_{n}}{1+pv_{n}+qu_{n}}$

Let T be the map on $\mathbb{C}\times\mathbb{C}$ to itself defined by

$T\left(\begin{array}[]{c}u\\\ v\\\
\end{array}\right)=\left(\begin{array}[]{c}v\\\ \frac{v}{1+pv+qu}\\\
\end{array}\right)$

Then $\left(\begin{array}[]{c}\phi\\\ \psi\\\ \end{array}\right)$ is a fixed
point of $T^{2}$, the second iterate of $T$.

$T^{2}\left(\begin{array}[]{c}u\\\ v\\\
\end{array}\right)=\left(\begin{array}[]{c}\frac{v}{1+pv+qu}\\\ \\\
\frac{\frac{v}{1+pv+qu}}{1+p\frac{v}{1+pv+qu}+qv}\\\ \end{array}\right)$

$T^{2}\left(\begin{array}[]{c}u\\\ v\\\
\end{array}\right)==\left(\begin{array}[]{c}g(u,v)\\\ h(u,v)\\\
\end{array}\right)$

where $g(u,v)=\frac{v}{1+pv+qu}$ and
$h(u,v)=\frac{\frac{v}{1+pv+qu}}{1+p\frac{v}{1+pv+qu}+qv}$. Clearly the two
cycle is locally asymptotically stable when the eigenvalues of the Jacobian
matrix $J_{T^{2}}$, evaluated at $\left(\begin{array}[]{c}\phi\\\ \psi\\\
\end{array}\right)$ lie inside the unit disk.

We have,

$J_{T^{2}}\left(\begin{array}[]{c}\phi\\\ \psi\\\
\end{array}\right)=\left(\begin{array}[]{cc}\frac{\delta g}{\delta
u}(\phi,\psi)&\frac{\delta g}{\delta v}(\phi,\psi)\\\ \\\ \frac{\delta
h}{\delta u}(\phi,\psi)&\frac{\delta h}{\delta v}(\phi,\psi)\\\
\end{array}\right)$

where $\frac{\delta g}{\delta
u}(\phi,\psi)=-\frac{q\psi}{(1+q\phi+p\psi)^{2}}$ and $\frac{\delta g}{\delta
v}(\phi,\psi)=\frac{1+q\phi}{(1+q\phi+p\psi)^{2}}$

$\frac{\delta h}{\delta
u}(\phi,\psi)=-\frac{q\psi(1+q\psi)}{\left(1+\psi+p\psi+q\left(\phi+\psi+q\phi\psi+p\psi^{2}\right)\right)^{2}}$
and $\frac{\delta h}{\delta
v}(\phi,\psi)=\frac{1+q\left(\phi-p\psi^{2}\right)}{\left(1+\psi+p\psi+q\left(\phi+\psi+q\phi\psi+p\psi^{2}\right)\right)^{2}}$

Now, set

$\chi=\frac{\delta g}{\delta u}(\phi,\psi)+\frac{\delta h}{\delta
v}(\phi,\psi)=-\frac{q\psi}{(1+q\phi+p\psi)^{2}}+\frac{1+q\left(\phi-p\psi^{2}\right)}{\left(1+\psi+p\psi+q\left(\phi+\psi+q\phi\psi+p\psi^{2}\right)\right)^{2}}$
$\lambda=\frac{q^{2}\psi^{2}}{(1+q\phi+p\psi)\left(1+\psi+p\psi+q\left(\phi+\psi+q\phi\psi+p\psi^{2}\right)\right)^{2}}$

In particular for the prime period $2$ solution,

$\left\\{\phi\to\frac{1}{0.5p-0.5q-0.5\sqrt{p^{2}-q^{2}}},\psi\to-\frac{1}{q}+\frac{\sqrt{p^{2}-q^{2}}}{pq-q^{2}}\right\\}$,
we shall see the local asymptotic stability for some example cases of
parameters $p$ and $q$. The general form of $\chi$ and $\lambda$ would very
complected. Consider the prime period two solution of the difference equation
(7), $\phi\to 0.365026+0.263198i,\psi\to-0.412345+0.131124i$ corresponding two
the parameters $p\to\frac{1}{5}+3i$ and $q\to\frac{3}{5}+5i$.

In this case, $\left\lvert\chi\right\rvert=0.0287948$ and
$\left\lvert\lambda\right\rvert=0.0000431717$. Therefore, by the Linear
Stability theorem ($|\chi|<1+|\lambda|<2$) the prime period $2$ solution
$\phi\to 0.365026+0.263198i,\psi\to-0.412345+0.131124i$ is _locally
asymptotically stable_.

###  5.3 Prime Period Two Solutions of Eq. (8)

Let $\ldots,\phi,~{}\psi,~{}\phi,~{}\psi,\ldots$, $\phi\neq\psi$ be a prime
period two solution of the difference equation
$w_{n+1}=\frac{w_{n-1}}{p+qw_{n}+w_{n-1}}$. Then
$\phi=\frac{\phi}{1+p\psi+q\phi}$ and $\psi=\frac{\psi}{1+p\phi+q\psi}$. This
two equations lead to the set of solutions (prime period two) except the
equilibriums as $\left\\{\phi\to 0,\psi\to 1-p\right\\}$.

Let $\ldots,\phi,~{}\psi,~{}\phi,~{}\psi,\ldots$, $\phi\neq\psi$ be a prime
period two solution of the equation (8). We set

$u_{n}=w_{n-1}$ $v_{n}=w_{n}$

Then the equivalent form of the difference equation (8) is

$u_{n+1}=v_{n}$ $v_{n+1}=\frac{u_{n}}{p+qv_{n}+u_{n}}$

Let T be the map on $\mathbb{C}\times\mathbb{C}$ to itself defined by

$T\left(\begin{array}[]{c}u\\\ v\\\
\end{array}\right)=\left(\begin{array}[]{c}v\\\ \frac{u}{p+qv+u}\\\
\end{array}\right)$

Then $\left(\begin{array}[]{c}\phi\\\ \psi\\\ \end{array}\right)$ is a fixed
point of $T^{2}$, the second iterate of $T$.

$T^{2}\left(\begin{array}[]{c}u\\\ v\\\
\end{array}\right)=\left(\begin{array}[]{c}\frac{u}{p+qv+u}\\\ \\\
\frac{\frac{u}{p+qv+u}}{p+q\frac{u}{p+qv+u}+v}\\\ \end{array}\right)$

$T^{2}\left(\begin{array}[]{c}u\\\ v\\\
\end{array}\right)==\left(\begin{array}[]{c}g(u,v)\\\ h(u,v)\\\
\end{array}\right)$

where $g(u,v)=\frac{u}{p+qv+u}$ and
$h(u,v)=\frac{\frac{u}{p+qv+u}}{p+q\frac{u}{p+qv+u}+v}$. Clearly the two cycle
is locally asymptotically stable when the eigenvalues of the Jacobian matrix
$J_{T^{2}}$, evaluated at $\left(\begin{array}[]{c}\phi\\\ \psi\\\
\end{array}\right)$ lie inside the unit disk.

We have,

$J_{T^{2}}\left(\begin{array}[]{c}\phi\\\ \psi\\\
\end{array}\right)=\left(\begin{array}[]{cc}\frac{\delta g}{\delta
u}(\phi,\psi)&\frac{\delta g}{\delta v}(\phi,\psi)\\\ \\\ \frac{\delta
h}{\delta u}(\phi,\psi)&\frac{\delta h}{\delta v}(\phi,\psi)\\\
\end{array}\right)$

where $\frac{\delta g}{\delta
u}(\phi,\psi)=\frac{p+q\psi}{(p+\phi+q\psi)^{2}}$ and $\frac{\delta g}{\delta
v}(\phi,\psi)=-\frac{q\phi}{(p+\phi+q\psi)^{2}}$

$\frac{\delta h}{\delta
u}(\phi,\psi)=\frac{(1+\psi)(p+q\psi)}{\left(p+\phi+q\phi+(p+q+\phi)\psi+q\psi^{2}\right)^{2}}$
and $\frac{\delta h}{\delta
v}(\phi,\psi)=-\frac{\phi(p+q+\phi+2q\psi)}{\left(p+\phi+q\phi+(p+q+\phi)\psi+q\psi^{2}\right)^{2}}$

Now, set

$\chi=\frac{\delta g}{\delta u}(\phi,\psi)+\frac{\delta h}{\delta
v}(\phi,\psi)=\frac{p+q\psi}{(p+\phi+q\psi)^{2}}-\frac{\phi(p+q+\phi+2q\psi)}{\left(p+\phi+q\phi+(p+q+\phi)\psi+q\psi^{2}\right)^{2}}$
$\lambda=-\frac{\phi(p+q\psi)}{(p+\phi+q\psi)\left(p+\phi+q\phi+(p+q+\phi)\psi+q\psi^{2}\right)^{2}}$

For the prime period $2$ solution, $\phi\to 0,\psi\to 1-p$,
$\chi=\frac{1}{p+(1-p)q}$ and $\left\lvert\lambda\right\rvert=0$. Therefore,
by the Linear Stability theorem ($|\chi|<1+|\lambda|<2$) the prime period $2$
solution is _locally asymptotically stable_ if and only if
$\left\lvert\frac{1}{p+(1-p)q}\right\rvert<1$. It turns out that $\left\lvert
q\right\rvert>\frac{1-\left\lvert p\right\rvert}{\left\lvert
1-p\right\rvert}$. In other words, the condition reduces to $\left\lvert
p\right\rvert<1$ and $\left\lvert q\right\rvert>1$ which is same condition as
it was for the real set up.

##  6 Chaotic Solutions

This is something which is absolutely new feature of the dynamics of the
difference equation (1) which did not arise in the real set up of the same
difference equation. Computationally we have encountered some chaotic solution
of the difference equation (8) for some parameter values which are given in
the following Table. 1.
The method of Lyapunov characteristic exponents serves as a useful tool to
quantify chaos. Specifically Lyapunav exponents measure the rates of
convergence or divergence of nearby trajectories. Negative Lyapunov exponents
indicate convergence, while positive Lyapunov exponents demonstrate divergence
and chaos. The magnitude of the Lyapunov exponent is an indicator of the time
scale on which chaotic behavior can be predicted or transients decay for the
positive and negative exponent cases respectively. In this present study, the
largest Lyapunov exponent is calculated for a given solution of finite length
numerically [10].
From computational evidence, it is arguable that for complex parameters $p$
and $q$ which are stated in the following table the solutions are chaotic for
every initial values.

Parameters $p$, $q$ |  Interval of Lyapunav exponent
---|---
$p=(0.2037,0.5444)$, $q=(0.8749,0.1210)$ |  $(0.3215,1.6235)$
$p=(0.4933,0.7018)$, $q=(0.8878,0.0551)$ |  $(1.062,2.021)$
$p=(0.7840,0.4867)$, $q=(0.4648,0.1313)$ |  $(0.6256,1.314)$
$p=(0.2308,0.6580)$, $q=(0.5629,0.2818)$ |  $(1.373,2.325)$
Table 1: Chaotic solutions of the equation (8) for different choice of
parameters and initial values.

The chaotic trajectory plots including corresponding complex plots are given
the following Fig. 1.

---

Figure 1: Chaotic Trajectories of the equation (8) of four different cases as
stated in Table 1.

In the Fig. 1, for each of the four cases ten different initial values are
taken and plotted in the left and in the right corresponding complex plots are
given. From the Fig. 1, it is evident that for the four different cases the
basin of the chaotic attractor is neighbourhood of the centre $(0,0)$ of
complex plane.

##  7 Some Interesting Nontrivial Problems

######  Open Problem 7.1.

Does the difference equation have higher order periodic cycle? If so, what is
the highest periodic cycle?

######  Open Problem 7.2.

Find out the set of all parameters $p$ and $q$ for which the difference
equation (8) has chaotic solutions.

######  Open Problem 7.3.

Find out the subset of the $\mathbb{D}$ of all possible initial values $z_{0}$
and $z_{1}$ for which the solutions of the difference equation are chaotic for
any complex parameters $p$ and $q$. Does the neighbourhood of $(0,0)$ is
global chaotic attractor? If not, are there any other chaotic attractors?

##  8 Future Endeavours

In continuation of the present work the study of the difference equation
${z_{n+1}=\frac{\alpha_{n}+\beta_{n}z_{n}+\gamma_{n}z_{n-1}}{A_{n}+B_{n}z_{n}+C_{n}z_{n-1}}}$
where $\alpha_{n}$, $\beta_{n}$, $\gamma_{n}$, $A_{n}$, $B_{n}$ and $C_{n}$
are all convergent sequence of complex numbers and converges to $\alpha$,
$\beta$, $\gamma$, $A$, $B$ and $C$ respectively is indeed would be very
interesting and that we would like to pursue further. Also the most
generalization of the present rational difference equation is

${z_{n+1}=\frac{\alpha+\beta z_{n-l}+\gamma z_{n-k}}{A+Bz_{n-l}+Cz_{n-k}}}$

where $l$ and $k$ are delay terms and it demands similar analysis which we
plan to pursue in near future.

## Acknowledgement

The author thanks _Dr. Pallab Basu_ for discussions and suggestions.

## References

  * [1] Saber N Elaydi, Henrique Oliveira, José Manuel Ferreira and João F Alves, _Discrete Dyanmics and Difference Equations, Proceedings of the Twelfth International Conference on Difference Equations and Applications_ , World Scientific Press, 2007.
  * [2] S. Atawna, R. Abu-Saris, E. S. Ismail, and I. Hashim, Stability of Nonhyperbolic Equilibrium Solution of Second Order Nonlinear Rational Difference Equation, _Journal of Difference Equations_ Volume 2015 (2015), Article ID 486985.
  * [3] Ch. G. Philos, I. K. Purnaras, and Y. G. Sficas, Global attractivity in a nonlinear difference equation, _Appl. Math. Comput._ 62(1994), 249-258.
  * [4] E. Camouzis and G. Ladas, _Dynamics of Third Order Rational Difference Equations; With Open Problems and Conjectures_ , Chapman & Hall/CRC Press, 2008.
  * [5] N. K. Govil and Q, I. Rahman, On the Eneström-Kakeya theorem, _Tohoku Mathematical Journal_ , 20(1968), 126-136.
  * [6] E. Camouzis, E. Chatterjee, G. Ladas and E. P. Quinn, On Third Order Rational Difference Equations – Open Problems and Conjectures, _Journal of Difference Equations and Applications_ , 10(2004), 1119 – 1127.
  * [7] G. T. Cargo and O. Shisha, Zeros of polynomials and fractional order differences of their coefficients, _Journ. Math. Anal. Appl._ , 7 (1963), 176-182.
  * [8] A. Joyal, G. Labelle and Q.I.Rahaman, On the location of zeros of polynomials, _Canad. Math. Bull._ , 10 (1967), 53-63.
  * [9] P. V. Krishnaiah, On Kakeya’s theorem, _Journ. London Math. Soc._ , 30 (1955), 314-319.
  * [10] A. Wolf, J. B. Swift, H. L. Swinney and J. A. Vastano, Determining Lyapunov exponents from a time series Physica D, 126(1985), 285-317.
  * [11] E.A. Grove, Y. Kostrov and S.W. Schultz, On Riccati Difference Equations With Complex Coefficients, _Proceedings of the 14th International Conference on Difference Equations and Applications, Istanbul, Turkey_ , March 2009.
  * [12] Sk. S. Hassan, E. Chatterjee, Dynamics of the equation $\displaystyle{z_{n+1}=\frac{\alpha+\beta z_{n}}{A+z_{n-1}}}$ in the Complex Plane, _Communicated to Computational Mathematics and Mathematical Physics, Springer_ , 2014
  * [13] M.R.S. Kulenovi$\acute{c}$ and G. Ladas, _Dynamics of Second Order Rational Difference Equations; With Open Problems and Conjectures_ , Chapman & Hall/CRC Press, 2001.
  * [14] V.L. Kocic and G. Ladas, _Global Behaviour of Nonlinear Difference Equations of Higher Order with Applications_ , Kluwer Academic Publishers, Dordrecht, Holland, 1993.
  * [15] Sk. S. Hassan, P. Basu, Complex Dynamics of the Difference Equation $z_{n+1}=\frac{\alpha}{z_{n}}+\frac{\beta}{z_{n-1}}$, arXiv:1502.06469 [math.DS].