Finite dynamical systems are dynamical systems on finite sets. Examples include cellular automata and Boolean networks, (e.g., [1]) with applications in many areas of science and engineering (e.g., [2, 3]), and more recently in computational biology (e.g., [4–6]). A common question in all of these applications is how to analyze the dynamics of the models without enumerating all state transitions. This paper presents partial solutions to this problem.

Because of technological advances such as DNA microarrays, it is possible to measure gene transcripts from a large number of genes. It is desirable to have efficient mathematical methods to extract information about regulatory iterations between genes from repeated measurements of gene transcript concentrations.

One piece of information about regulatory iterations of interest is when the dynamics reaches a steady state. In the words of Fuller (see [7]): “this paradigm closely parallels the goal of professionals who aim to understand the flow of molecular events during the progression of an illness and to predict how the disease will develop and how the patient will respond to certain therapies.”

The work of Fuller et al. [7] serves as an example. When the gene expression profile of human brain tumors was analyzed, these were divided into three classes—high grade, medium grade, and low grade. A key gene expression event was identified, which was a high expression of insulin-like growth factor binding protein 2 (IGFBP2) occurring only in high-grade brain tumors. It can be assumed that gene expression events were initiated at some stages in low-level tumors and may have led to the state when IGFBP2 is activated. The activation of IGFBP2 can be understood to be a steady state. If we model the kinetics and construct a model that reconstructs the genetic regulatory network that activates during the brain tumor process, then we may be able to predict the convergence of events that lead to the activation of IGFBP2. In the same way, we also want to know what happens in the next step following the activation of IGFBP2. Our goal is to develop tools that will enable this type of analysis in the case of modeling gene regulatory networks by means of discrete dynamical systems.

The use of polynomial dynamical systems to model biological phenomena, in particular gene regulatory networks, have proved to be as valid as continuous models. Laubenbacher and Stigler (see [6]) point out, for example, that most ordinary differential equations models cannot be solved analytically and that numerical solutions of such time-continuous systems necessitate approximations by time-discrete systems, so that ultimately, the two types of models are not that different.

Once a gene regulatory network is modeled, in our case by finite fields, or by finitely generated modules, we obtain a finite dynamical system. Our goal is to determine if the dynamical system represents a steady-state gene regulatory network (i.e., if every state eventually enters a steady state). This is a crucial task. Shmulevich et al. (see [8]) have shown that the steady-state distribution is necessary in order to compute the long term influence that is a measure of gene impact over other genes.

The rest of the paper is organized as follows. In Section 2 we give some basic definitions and facts about finite dynamical systems and their associated state spaces. In Section 3 we discuss multivariate and univariate finite field models for genetic networks and show that they are equivalent. Each of the models can be converted to the other by a discrete Fourier transform. Section 4 is devoted to fixed point systems. We give a new definition of linear finite dynamical systems and give necessary and sufficient conditions for such a system to be a fixed point system. We review results concerning monomial fixed point systems and show how our results concerning linear systems can be used to determine when a monomial finite dynamical system over an arbitrary finite field is a fixed point system. We show how fixed points can be determined in the univariable model by solving a polynomial equation over a finite field and we give a necessary condition for a finite dynamical system to be a fixed point system. Finally, in Section 5 we discuss some implementation issues.

A finite dynamical system (fds) is an ordered pair (X, f) where X is a finite set and f is a function that maps X into itself, that is, f : X → X. The state space of an fds (X, f) is a digraph (i.e., directed graph) whose nodes are labeled by the elements of X and whose edges consist of all ordered pairs (x, y) X × X such that f(x) = y. We say that two finite dynamical systems are isomorphic if there exists a graph isomorphism between their state spaces.

Let G = (V, E) be a digraph. A path in G of the form (v₁, v₂), (v₂, v₃), …, (v_{n − 1}, v_n), (v_n, v₁), where v₁, v₂, …, v_n are distinct members of V is a cycle of length n. We define a tree to be a digraph T = (V, E) which has a unique node v₀, called the root of T, such that (a) (v₀, v₀) E, (b) for any node v ≠ v₀, there is a path from V to v₀, (c) T has no “semicycles” (i.e., alternate sequence of nodes and edges v₁, x₁, v₂, …, x_n, v_n+1, n ≠ 0, where v₁ = v_n+1 and each x_i is (v_i, v_i+1 or v_i+1, v_i)) other than the trivial one (v₀, (v₀, v₀), v₀). (Such a tree with the edge (v₀, v₀) deleted is sometimes called an “in-tree” with “sink” v₀.)

Let T be a tree, let nT = _{i = 1}ⁿT_i be the union of n copies T₁, T₂, …, T_n of T, and let r_i be the root of T_i. Define T⁽ⁿ⁾ to be the digraph obtained from nT by deleting the edges (r_i, r_i), i = 1, 2, …, n, and adjoining the edges (r_i, r_j), i, j = 1, 2, …, n, where j = i + 1 mod n. We call T⁽ⁿ⁾ an n-cycled tree. Note that by definition, every tree is an n-cycled tree with n = 1. Note also that by definition, a digraph T_r = ({r}, {r, r}) consisting of a single trivial cycle is a tree and hence every cycle of length n is isomorphic to nT_r and hence is an n-cycled tree.

The product of two digraphs G₁ = (V₁, E₁) and G₂ = (V₂, E₂), denoted G₁× G₂, is the digraph G = (V, E) where V = V₁× V₂ (the Cartesian product of V₁ by V₂) and E = {((x₁, y₁), (x₂, y₂)) V × V : (x₁, x₂) E₁ and (y₁, y₂) E₂}. The following facts follow easily from the definitions. Lemmas 1, 3, and 4 have been noted in [9].

The state space of an fds is the disjoint union of cycled trees.

Of special interest are those fds whose state space consists entirely of trees. Such an fds is called a fixed point system (fps).

For any finite set X we call f : X → X nilpotent if there exists a unique x₀ X such that f^k(X) = x₀ for some positive integer k.

The state space of an fds (X, f) is a tree if and only if f is nilpotent. Hence (X, f) is an fps if f is nilpotent.

Suppose that the state space (X, f) is a tree with root x₀ and height k. Then f^k(x) = x₀ for all x X and x₀ is the only node with this property. Hence f is nilpotent. Conversely, if f is nilpotent and f^k(X) = x₀, then by Lemma 1, the state space consists of an n-cycled tree and since x₀ is unique, n = 1.

Consider the fds (F₂³, f), where f : F³₂ → F³₂ is defined by f(x, y, z) = (y, 0, x) and F₂ is the binary field. In this case f is a nilpotent function. The state space of (F₂³, f) is a tree whose state space is shown in Figure 1.

Figure 1 State space of (F23, f), where f(x, y, z) = (y, 0, x) over F2.

The state space of an fds (X, f) is the union of cycles if and only if f is one-to-one.

The product of a tree and a cycle of length l is a cycled tree whose cycle has length l.

A finite dynamical system constitutes a very natural discrete model for regulatory processes (see [10]), in particular genetic networks. Experimental data can be discretized into a finite set X of expression levels. A network consisting of n genes is then represented by an fds (Xⁿ, f). The dynamics of the network is described by a discrete time series

Special cases of the finite dynamical model are the Boolean model and finite field models. In the Boolean model, either a gene can affect another gene or not. In a finite field model, one is able to capture graded differences in gene expression. A finite field model can be considered as a generalization of the Boolean model since each Boolean operation can be expressed in terms of the sum and product in Z₂. In particular,

.

(2)

Two types of finite field models have emerged, the multivariate model [6] and the univariable model [11]. The multivariate model is given by the fds (F_qⁿ, f), where F_qⁿ represents the set of n-tuples over the finite field F_q with q elements. Each coordinate function f_i gives the next state of gene i, given the states of the other genes. The univariate model is given by the fds (F_qⁿ, f). In this case, each value of f represents the next states of the n genes, given the present states.

The two types of finite field models can be considered equivalent in the following sense.

An fds (X, f) is equivalent to an fds (Y, g) if there is an epimorphism : X → Y such that f = g .

It is easy to see that if two fds's are equivalent, then their state spaces are the same up to isomorphism. We can show that for any n-dimensional dynamical system (Fⁿ_q, f) there is an equivalent one-dimensional system (F_qⁿ, g). To see this, consider a primitive element α of F_qⁿ, that is, a generator of the multiplicative group of F_qⁿ − {0}. Then there is a natural correspondence between Fⁿ_q and F_qⁿ, given by

α(x0,…, xn − 1) = x0 + x1α + x2α2 + + xn − 1αn − 1.

(3)

a F_qⁿ

y_i F_q

a = y₀ + y₁α + y₂α² + + y_{n − 1}α^{n − 1}

g _α = _α f

One important consideration in choosing an appropriate finite field model for a genetic network is the complexity of the needed computational tasks. For example, the evaluation of a polynomial in n variables over F_q, q prime, can be done with O(qⁿ/n) operations (see [12]) and hence, evaluating f in all n of its coordinates is O(qⁿ), the same number of operations needed for the evaluation of a univariate polynomial over F_qⁿ. However, the complexity of the comparison of two values in F_qⁿ is O(n), whereas the complexity of the comparison of two values in F_qⁿ, represented as described below, is O(1).

Arithmetic in F_qⁿ, q prime, is integer arithmetic modulo q. Arithmetic in F_qⁿ is efficiently performed by table lookup methods, as shown below. Nonzero elements of F_qⁿ are represented by powers of a primitive element α. Multiplication is then performed by adding exponents modulo q^{n − 1}. For addition we make use of a precomputed table of values defined as follows. Every nonozero element of F_qⁿ has a unique representation in the form 1 + αⁱ and the unique number z(i), 0 ≤ z(i) ≤ q^{n − 2}, such that 1 + αⁱ = α^z(i) is called the Zech log of i. Note that for a ≤ b, α^a + α^b = α^a(1 + α^{b − a}) = α^{a + z(b − a)modqn − 1}. Addition is thus performed by adding one exponent to the Zech log of the difference, which is found in a precomputed table. In order to construct a table of Zech logs for F_qⁿ, we first need a primitive polynomial, which can be found in any one of various tables (e.g., [13]).

Let us construct a table of Zech logs for F₃₂ using the primitive polynomial x⁵ + x² + 1. Thus, we have α⁵ = α² + 1, where α is a root of x⁵ + x² + 1. Continuing to compute the powers and making use of this fact, we have α⁶ = α³ + α, α⁷ = α⁴ + α², α⁸ = α⁵ + α³ = α³ + α² + 1,…, α³¹ = 1. Now use these results to compute for each i = 1,…, 30, the number z(i) such that αⁱ + 1 = α^z(i). For example, since α⁵ = α² + 1, we have α⁵ + 1 = α² and so z(5) = 2, and so forth. See Table 1.

Table 1 Zech Logs for F32.

Usually it is most convenient to choose the most appropriate model at the outset. However, at the cost of computing all possible values of the map, it is possible to convert one model to the other. The rest of this section is devoted to developing such an algorithm.

Let F be a field and let α F be an element of order d, that is, α^d = 1 and no smaller power of α equals 1. The discrete Fourier transform (DFT) of blocklength d over F is defined by the matrix T = [α^ij], i, j = 0, 1,…, d − 1. The inverse discrete Fourier transform is given by T ^{− 1} = d ^{− 1}[α ^{− ij}], i, j = 0, 1,…, d − 1, where d ^{− 1} denotes the inverse of the field element d = 1 + 1 + + 1 (d times).

It is easy to show that TT ^{− 1} = I_d, where I_d denotes the d × d identity matrix (see, e.g., [14]). Now an element in F_q, is of order d if and only if d divides q − 1. Thus, for every finite field F_q there is a DFT over F_q with block length q − 1 which is defined by [α^ij], i, j = 0, 1,…, q − 2, where α is a primitive element of F_q. We denote such a DFT by T_{q, α}.

Let B₀ = (_α f)(0,…, 0) and for each i = 1, 2,…, qⁿ − 1, let B_i = (_α f) (a_{0, i}, a_{1, i}, …, a_{n − 1, i}) where α is a primitive element of F_qⁿ and where a_{n − 1, i}α^{n − 1} + + a_{1, i}α + α_{0, i} = α^{i − 1}. Then g is given by the polynomial

Aqn − 1xqn − 1 + Aqn − 2xqn − 2 + + A1x + A0,

(4)

(5)

For each i = 0, 1,…, qⁿ − 2, we have

.

(6)

g(αi) = Aqn − 1(αi)qn − 1 + Aqn − 2(αi)qn − 2 + + A1αi + A0, i = 0, 1,…, qn − 2.

(7)

Bi+1 − A0 = g(αi) = Aqn − 1(αi)qn − 1 + Aqn − 2(αi)qn − 2 + + A1αi, i = 0, 1,…, qn − 2.

(8)

,

(9)

We illustrate the algorithm given by Theorem 1 with an example.

A recent application involves the study and creation of a model for lac operon [15]. When the bacteria E. Coli is in an environment with lactose, then the lac operon turns on the enzymes that are needed in order to degrade lactose. These enzymes are beta-galactisidase, Lactose Permease, and Thiogalactoside transectylase. In [15], a continuous model is proposed that measures the rate of change in the concentration of these enzymes as well as the concentration of mRNA and intracellular lactose. In [16, 17], Laubenbacher and Stigler provide a discrete model for the lac operon given by

In our case, we choose α to be a zero of x⁵ + x² + 2. Next, we compute B_i, i = 0, 1,…, 31. By definition B₀ = _α(f(0, 0, 0, 0, 0)) = _α(0, 0, 0, 0, 0) = 0 and B_i = _α(f(a_{0, i}, a_{1, i}, a_{2, i}, a_{3, i}, a_{4, i})) where α^{i − 1} = a_{0, i} + a_{1, i}α + a_{2, i}α² + a_{3, i}α³ + a_{4, i}α⁴ for i = 1, 2,…, 31.

So, for example,

.

(11)

Multiplying by the 31 × 31 matrix T₃₂, α = [α^ij], 0 ≤ i, j ≤ 30, we obtain [A₃₁, A₃₀,…, A₁] and hence the equivalent univariate polynomial, which is g(x) = x + α²²x² + α²x³ + α¹¹x⁴ + α²⁰x⁵ + x⁶ + α¹²x⁷ + α²⁵x⁸ + α²⁰x⁹ + α²⁰x¹⁰ + α²x¹¹ + α⁵x¹² + α⁵x¹³ + α²³x¹⁴ + α⁷x¹⁶ + α²⁷x¹⁷ + α²⁰x¹⁸ + α⁶x¹⁹ + α²⁰x²⁰ + α²⁷x²¹ + x²² + α²⁷x²⁴ + x²⁵ + α²⁷x²⁶ + α¹⁶x²⁸.

As previously mentioned, the complexity of evaluating a polynomial in n variables over a finite field F_q is O(qⁿ/n). The complexity of evaluating f in all of its n coordinates is thus O(qⁿ) and the complexity of evaluating f in all points of F_qⁿ is thus O(q²ⁿ). The computation of the matrix-vector product in (2) involves O(q²ⁿ) operations over the field F_qⁿ. However, using any one of the number of classical fast algorithms, such as Cooley-Tukey (see, e.g., [14]), the number of operations can be reduced to O(qⁿn).

A fixed point system (fps) is defined to be an fds whose state space consists of trees, that is, contains no cycles other than trivial ones (of length one). The fixed point system problem is the problem of determining of an fds whether or not it is an fps. Of course one such method would be the brute force method whereby we examine sequences determined by successive applications of the state map to determine if any such sequence contains a cycle of length greater than one. The worst case occurs when the state space consists of one cycle. Consider such a multivariate fds (F_qⁿ, f). In order to recognize a maximal cycle, f(a₁, a₂,…, a_n), f²(a₁, a₂,…, a_n),…, such an approach would require backtracking at each step in order to compare the most recent value fⁱ(a₁, a₂,…, a_n) with all previous values. An evaluation requires O(qⁿ) operations, there are qⁿ such evaluations, and a comparison of two values requires n steps. The complexity of the complete process is thus O(q²ⁿ+n²) = O(q²ⁿ).

To date, all results concerning the fixed point system problem are characterized in terms of multivariate fds. A solution to the fixed point system problem consists of characterizing such an fds (F_qⁿ, f) in terms of the structure of f. Ideally, such conditions should be amenable to implementations in polynomial time in n. In a recent work, Just [18] claims that if the class of regulatory functions contains the quadratic monotone functions x_i x_j and x_i ∧ x_j, then the fixed point problem for Boolean dynamical systems is NP-hard. In view of this result, it is unlikely that we can achieve the goal of a polynomial time solution to the fixed point problem, at least in the general case. However, the question arises if the above O(q²ⁿ) result can be improved (to say O(qⁿ)) and also what are special cases of the fixed point problem that have polynomial time solutions.

In this section we give a polynomial solution to the special case of the fixed point problem for linear finite dynamical systems, we review known results for the nonlinear case, and we point out how our results concerning the general linear case for fds over finitely generated modules give a more complete solution to the case of monomial finite field dynamical systems over arbitrary finite fields. We conclude by proposing a new approach to the problem via univariate systems, we give an algorithm for determining the fixed points of a univariate system, and we give a necessary condition for a univariate fds to be an fps.

In the multivariate case, a linear system over a finite field is represented by an fds (F_qⁿ, f) where f can be represented by an n × n matrix A over F_q. The fixed points of a multivariate fds (F_q^m, A) are simply the solutions to the homogeneous system of equations (A − I)x = 0.

In the finite field model for genetic networks, we assume that the number of states of each gene is a power of a prime. However, we will give a more general model that eliminates this assumption.

A module M over a ring R is said to be finitely generated if there exists a set of elements {s₁, s₂,…, s_n} M such that M = {r₁s₁ + r₂s₂ + + r_ns_n r_i R}. Finitely generated modules are generalized vector spaces. Examples are F_qⁿ and the set Z_mⁿ of n tuples over the ring of integers modulo an arbitrary integer m.

A linear finite dynamical module system (lfdms) consists of an ordered pair (M(R), f) where M(R) is a finitely generated module over a finite commutative ring R with unity and f : M(R) → M(R) is linear. Let (M₁(R), f₁) and (M₂(R), f₂) be lfdms. We define the direct sum of (M₁(R), f₁) and (M₂(R), f₂) to be the fds (M₁ M₂, f₁ f₂) where M₁ M₂ is the direct sum of the modules M₁(R) and M₂(R) and f₁ f₂ : M₁ M₂ → M₁ M₂ is defined by (f₁ f₂)(u + v) = f₁(u) + f₂(v), for each u M₁(R) and v M₂(R). The state space of the direct sum is related to the component fds as follows.

Let G₁ be the state space of the lfdms (M₁(R), f₁) and let G₂ be the state space of the lfdms (M₂(R), f₂). Then the state space of the direct sum of (M₁(R), f₁) and (M₂(R), f₂) is G₁ × G₂.

This result has been noted in [9] for lfdms over fields.

We use the following well-known result (see, e.g., [19]) in order to establish necessary and sufficient conditions for an lfdms to be a fixed point system.

Let (M(R), f) be an lfdms. Then there exist an integer n > 0 and submodules N and P satisfying

• N = fⁿ(M(R)),
• P = f⁻ⁿ(0),
• (M(R), f) = (N(R), f₁) (P(R), f₂), f₁ = f|_N (the restriction of f to N) is invertible and f₂ = f|_P is nilpotent.

Let (M(R), f) be an lfdms and let N be defined as above. Then (M(R), f) is a fixed point system if and only if either f is nilpotent or f|_N is the identity map.

By Fitting's lemma, we have (M(R), f) = (N(R), f|_N) (P(R), f|_P) where N = fⁿ(M(R)) and P = f⁻ⁿ(0). Suppose that f is nilpotent. Then by Lemma 2, the state space of (M(R), f) is a tree. Next suppose that f|_N is the identity. Then by Lemma 3, the state space of (M(N), f|_N) is a union of cycles each of length one and by Lemma 2, the state space of (M(P), f|_P) is a tree. Hence by Lemma 4, the state space of (M(R), f) is a union of trees and so (M(R), f) is a fixed point system.

Conversely, suppose that (M(R), f) is a fixed point system. Then the state space of (M(R), f) is a union U of trees. If U consists of only one tree, then by Lemma 3, f is nilpotent. Now suppose that U is the union of at least two trees. Since f is invertible on N, it is also one-to-one on N. By Lemma 2, the state space of (N(R), f|_N) is a union of cycles. Each of these cycles must be of length one. For if not, the state space of (M(R), f) would contain at least one n-cycled tree where n > 1, contradicting that (M(R), f) is a fixed point system.

Theorem 2 can be used to prove the following result, which is suggested in [20, 21].

A linear finite dynamical system (F_qⁿ, f) over a field is a fixed point system if and only if the characteristic polynomial of f is of the form x^n₀(x − 1)^n₁ and the minimal polynomial is of the form x^n₀′(x − 1)^n₁′ where n₁′ is either zero or one.

Suppose (F_qⁿ, f) is an fps. Then either f is nilpotent or f|_N is the identity. If f is nilpotent, then the characteristic polynomial of f is of the form x^n₀ and the minimal polynomial of f is of the form x^n₀′. If f|_N is the identity, then the characteristic polynomial of f|_N is of the form (x − 1)^n₁ and the minimal polynomial of f|_N is of the form (x − 1)^n₁′ where 0 ≤ q n₁′ ≤ q n₁. Furthermore, n₁′ ≤ q 1 since otherwise (F_qⁿ, f) would not be an fps [2]. Therefore the characteristic and minimal polynomials of f are of the desired forms.

Conversely, suppose that the characteristic polynomial of f is of the form x^n₀(x − 1)^n₁ and its minimal polynomial is of the form x^n₀′(x − 1)^n₁′, where n₀′ ≤ q n₀ and n₁′ is either zero or one. If n₀′ = 0, then the characteristic polynomial of f is (x − 1)^n₁ and so the minimal polynomial of f is x − 1, which implies that f is the identity and hence (F_qⁿ, f) is an fps. Next, suppose that n₀′ > 0. Then either n₁′ = 0 or n₁′ = 1. If n₁′ = 0, then the state space of (F_qⁿ, f) is a tree. If n₁′ = 1, then the state space of (F_qⁿ, f) is the product of a tree and cycles of length one and, hence the union of trees.

The corollary gives us a polynomial time algorithm to determine of a linear fds (Fⁿ_q, f), where f is given by an n × n matrix, whether or not it is an fps. The characteristic polynomial of f can be determined in time O(n³) using the definition. The minimal polynomial of f can be determined in time O(n³) using an algorithm of Storjohann [22]. Both polynomials can be factored in subquadratic time using an algorithm of Kaltofen and Shoup [23].

If f = (f₁, f₂,…, f_n) where each f_j is of the form x₁ ^ε_1jx₂ ^ε_2j x_n ^ε_nj, j = 1, 2,…, n, where each ε_ij belongs to the ring Z_q−1 of integers modulo q − 1, then (F_qⁿ, f) is called a monomial finite dynamical system. The log map of (F_qⁿ, f) is defined by the n × n matrix L_f = [ ε_ij], where 1 ≤ i, j ≤ n. The support map is defined by S_f = (h₁, h₂,…, h_n) where each h_i = x₁^δ_i1x₂^δ_i2 x_n^δ_in and where δ_ij is one if ε > 0 and is zero otherwise.

The following theorem was published in [25].

A monomial fds (F_qⁿ, f) is an fps if and only if (Z_q−1ⁿ, L_f) and (Z₂ⁿ, S_f) are fixed point systems.

Consider the monomial fds (F₅², f) where f = (xy, y) = (x¹y¹, x⁰y¹). The matrix L_f = over Z₄ is nonsingular and hence not nilpotent. Furthermore, the n of Theorem 2 is 1, N = Z₄, and L_f is not the identity. By Theorem 2, (Z₄, L_f) is not an fps and by the previous theorem, (F₅², f) is not an fps.

The problem of determining in polynomial time (in n) when an lfdms (Rⁿ, f) is an fps, where R is a finite ring, is open.

(F_qⁿ, g) has fixed points if and only if h(x) = gcd(g(x) − x, x^qn − x) ≠ 1 and in such a case, the fixed points are the zeros of h(x).

An element a of F_qⁿ is a fixed point of (F_qⁿ, g) if and only if a is a zero of g(x) − x. Since x^qn = x for all x F_qⁿ, x^qn − x contains all linear factors of the form x − a, a F_qⁿ and so a is a zero of g(x) − x if and only if x − a is a factor of both g(x) − x and x^qn − x, that is, if and only if it is a factor of h(x).

Lemma 7 gives us algorithms for determining whether or not a given univariate fds has fixed points and if so, a method to find all such points. For the first part, we note that the greatest common divisor of two univariate polynomials of degree no more than d can be determined using no more than (d log² d) operations [26]. Since g has degree at most qⁿ − 1, this means that the complexity for calculating h(x), that is, for determining whether or not a given univariate fds (F_qⁿ, g) is an fps O(n²qⁿ).

When h(x) ≠ 1, h(x) can be factored in order to determine the set of all fixed points. At worst, using the algorithm in [23], the complexity of determining the factors of h(x) is O(d^1.815n), where d is the degree of h(x). Clearly, d is less than or equal to the degree of g(x)which in practice is determined by experimental data (e.g., from microarrays) and thus considerably less than the total number of possible points qⁿ. If we assume that the degree of g(x) is not more than the square root of qⁿ, then d^1.815n ≤ n²qⁿ and the total complexity of the algorithm for determining all fixed points is thus O(n²qⁿ).

In contrast, the only known method for determining the fixed points of a multivariate fds (F_qⁿ, f) is the brute force method of enumerating all state transitions and for each value f(a₁, a₂,…, a_n) so generated, check to see if f(a₁, a₂,…, a_n) = (a₁, a₂,…, a_n). The number of operations in this method is O(q²ⁿ).

In many cases, the degree of h(x) of Lemma 7 is small and its zeros can be found by inspection or by only several trials and errors. The lac operon example illustrates this.

Let (F₃₂, g) be the fds describing the lac operon (Example 3). We have h(x) = gcd(g(x), x³² − x) = x⁴ + α²⁶x³ + α¹⁸x² = x²(x − α³)(x − α¹⁵) and thus the fixed points are x = 0, x = α³, and x = α¹⁵.

Lemma 7 also gives a necessary condition for an fds to be an fps, which for emphasis we state as a theorem.

With the notation of Lemma 7, if (F_q, g) is an fps, then h(x) ≠ 1.

If h(x) = 1, then (F_qⁿ, g) has no fixed points and all cycles are nontrivial. Hence by Lemma 7, (F_qⁿ, g) is not an fps.

The converse of Theorem 4 is not true.

Consider the fds (F₅, g) where g(x) = x³. Then h(x) = gcd(x⁵ − x, x³ − x) = x³ − x ≠ = 1, but (F₅, g) is not an fps. (see Figure 2).

Figure 2 State space of (F5, g), where g(x) = x3 over F5.

One of the difficulties of implementing algorithms for the multivariate model is the choice of data structures, which can, in fact, affect complexity. For example, no algorithm is known for factoring multivariate polynomials that runs in time polynomial in the length of the “sparse” representation. However, such an algorithm exists for the “black box” representation (see, e.g., [27]).

On the other hand, data structures needed for algorithms for the univariate model are well known and simple to implement. In this case, one can also take advantage of well-known methods used in cryptography and coding theory. Table lookup methods for carrying out finite field arithmetic are an example. By using lookup tables we can make arithmetic operations at almost no cost. However, for very large fields, memory space becomes a limitation. Ferrer [28] has implemented table lookup arithmetic for fields of characteristic 2 on a Hewlett-Packard Itanium machine with two 900MHz ia64 CPU modules and 4GB of RAM. On this machine, we can create lookup tables of up to 2²⁹ elements.

Multiplication is by far the most costly finite field operation and also the most often used, since other operations such as computing powers and computing inverses make use of multiplication. In other experiments on the Hewlett-Packard Itanium, Ferrer [28] takes advantage of machine hardware in order to implement a “direct” multiplication algorithm for F_2ⁿ that runs in time linear in n for n = 2 up to n = 63 [28]. Here the field size is limited by the word-length of the computer architecture.

For larger fields, we can make use of “composite” fields (see, e.g., [29]), that is, fields F_n where n is composite, say n = rs. Making use of the isomorphism of F_2^rs and F_(2^r)^s, we can use table lookup for a suitable “ground field” F_r and the direct method mentioned above for multiplication in the extension field F_(2^r)^s. Using the ground field F_2⁵ and selected values of s, Ferrer [28] obtains running time O(s²).

Still another approach to implement finite field arithmetic, that is, especially efficient for fields of characteristic 2, is the use of reconfigurable hardware or “field programmable gate arrays” (FPGAs). In [30], Ferrer, Moreno and the first author obtain a multiplication algorithm which outperforms all other known FPGA multiplication algorithms for fields of characteristic 2.

One piece of information that is of utmost interest when modeling biological events, in particular gene regulation networks, is when the dynamics reaches a steady state. If the modeling of such networks is done by discrete finite dynamical systems, such information is given by the fixed points of the underlying system. We have shown that we can choose between a multivariate and a univariate polynomial representation. Here we introduce a new tool, the discrete Fourier transform that helps us change from one representation to the other, without altering the dynamics of the system.

We provide a criterion to determine when a linear finite dynamical system over an arbitrary finitely generated module over a commutative ring with unity is a fixed point system. When a gene regulation network is modeled by a linear finite dynamical system we can then decide if such an event reaches a steady state using our results. When the finitely generated module is a finite field we can decide in polynomial time.

Gene regulation networks, as suggested in the literature, seem to obey very complex mechanisms whose rules appear to be of a nonlinear nature (see [31]). In this regard, we have made explicit some useful facts concerning fixed points and fixed point systems. We have given algorithms for determining when a univariate fds has at least one fixed point and how to find them. We have also given a necessary condition for a univariable fds to be a fixed point system. However, there are still much to be done and a number of open problems remain. In particular, what families of fds admit polynomial time algorithms for determining whether or not a given fds is an fps? This work is a first step towards the aim of designing theories and practical tools to tackle the general problem of fixed points in finite dynamical systems.

This work was partially supported by Grant number S06-GM08102 NIH-MBRS (SCORE). The figures in this paper were created using the Discrete Visualizer of Dynamics software, from the Virginia Bioinformatics Institute (http://dvd.vbi.vt.edu/network_visualizer). The authors are grateful to Dr. Oscar Moreno for sharing his ideas on the univariate model and composite fields.