A Continuous-Time Network Evolution Model Describing 2- and 3-Interactions

: A continuous-time network evolution model is considered. The evolution of the network is based on 2- and 3-interactions. 2-interactions are described by edges, and 3-interactions are described by triangles. The evolution of the edges and triangles is governed by a multi-type continuous-time branching process. The limiting behaviour of the network is studied by mathematical methods. We prove that the number of triangles and edges have the same magnitude on the event of non-extinction, and it is e α t , where α is the Malthusian parameter. The probability of the extinction and the degree process of a ﬁxed vertex are also studied. The results are illustrated by simulations. set { r 1 , p 1 , p 2 , b , c } , we present the conﬁdence intervals calculated from the number of edges being born ( E ) resp. being present ( ˜ E ) and from the number of triangles being born ( T ) resp. being present ( ˜ T ) up to time t = 14. The conﬁdence intervals containing the numerical Malthusian parameter ˆ α are highlighted with the ∗ symbol. We see that any conﬁdence interval is narrow, and it either contains ˆ α , or ˆ α is very close to the interval. These results show that the approximation is good for moderate values of t .


Introduction
Network theory has a vast literature. In the book of Barabási [1], the general aspects can be found, while the book of van der Hofstad [2] is devoted to the mathematical models. Any network can be considered as a graph. The nodes of the network are the vertices, and the connections are the edges of the graph. A most famous model is the preferential attachment model proposed by Albert and Barabási [3]. It is a discrete time network evolution model and it describes connections of two nodes. In real life, the meaning of connection can be any interaction or any cooperation.
There are models for cooperation more than two units. For example, Backhausz and Móri studied three-interactions in [4]. Their model is generalised for N-interactions by Fazekas and Porvázsnyik in [5]. Both of these papers consider cliques where, inside a team, all members cooperate. In some sense, the opposite of the cliques, i.e., star-like connections were considered by Fazekas and Perecsényi [6]. In [6], there is no cooperation between two peripheral members of the team but all of them cooperate with the central member of the team. Despite [3], in papers [4][5][6] the preferential attachment rule is used for certain subgraphs and not for vertices.
We mention that in [7] the Erdős-Rényi graph, the configuration model and the preferential attachment graph were studied when the population was split into two types. The mathematical tool of the analysis in [7] is the theory of multi-type branching processes.
There are several continuous-time network evolution models. Here, we list only some papers using continuous time branching processes. Early works in this direction are [8,9]. Recently, in [10], multi-type preferential attachment trees were studied. In [10], the results of [11] on multi-type continuous time branching processes were applied to describe the evolution of the network.
In this paper, we study a new network evolution model. The structure and the rules of the evolution of our model were inspired both by some everyday experiences and deep scientific results on motifs. On the one hand, we had in our mind activities and the Malthusian parameter. To prove Theorem 2, we use the underlying branching process counted with certain random characteristics and apply the asymptotic theorems of [11].
In Section 5, the generating functions are calculated. Using the generating functions, the probability of extinction are studied. In Section 6, the asymptotic behaviour of the degree of a fixed vertex is considered. Here, we again apply the asymptotic theorems of [11] but with other characteristics than in Section 4. In Section 7, we present some simulation results supporting our theorems. Our figures and tables show that the values obtained by simulation fit well to the theoretical results.
The proofs are based on known general results of multi-type continuous-time branching processes. Therefore, for the reader's convenience, in Section 8, we list several results on multi-type Crump-Mode-Jagers processes.
We mention that our model was presented in our conference paper [18]. In that paper some preliminary theoretical results were announced together with some numerical evidence but without mathematical proofs.

The Model
We study the following network evolution model. At the initial time t = 0 the network consists of one single object, this object can be either an edge or a triangle. This object is called the ancestor. During the evolution, this ancestor object produces offspring objects, which can be either edges or triangles. Then, these offspring objects produce their offspring objects and so on. The reproduction times of any fixed object, including the ancestor, are the occurrences in its own Poisson process with rate 1.
From the theory of branching processes, we apply the following usual assumptions. That is we suppose that the reproduction processes of different objects are independent. Moreover, we assume that the reproduction processes of the edges are independent copies of the reproduction process of the generic edge. Similarly, the reproduction processes of the triangles are independent copies of the reproduction process of the generic triangle.
First, we explain the evolution of the generic edge. A Poisson process Π 2 (t) with parameter 1 gives its reproduction times. At any jumping time of this Poisson process, a new vertex appears and it is connected to the generic edge with one or two edges. The probability that this new vertex is connected to the generic edge by one new edge is r 1 , where 0 ≤ r 1 ≤ 1. The other end point of this new edge is chosen from the two vertices of the generic edge uniformly at random. We see that in this case the generic edge produces always one new edge. The other case is that when the new vertex is connected to both vertices of the generic edge. Its probability is r 2 = 1 − r 1 . In this second case the offspring of the generic edge is a triangle consisting of the generic edge and the two new edges. We emphasize that in this last case the generic edge itself and the new triangle will produce offspring, but the two new edges are not substantive parts of the reproduction process, so they alone will not produce offspring.
The reproduction process of the generic triangle is similar. The Poisson process with rate 1 corresponding to the generic triangle is denoted by Π 3 (t), t ≥ 0. The jumping times of Π 3 (t) are the birth times of the generic triangle. At every birth time a new vertex is born and it joins to the existing graph so that it is connected to our generic triangle with 1, 2 or 3 edges. Denote by p j (j = 1, 2, 3) the probability that the new vertex is connected to j vertices of our generic triangle. The vertices of the generic triangle to be connected to the new vertex are chosen uniformly at random.
By the above definition of the evolution process, at each birth step we add precisely 1 new vertex. When the new vertex is connected to one vertex of the generic triangle, the generic triangle gives birth to one new edge. This event has probability p 1 . However, in the remaining two cases we count only the new triangles and not the new edges. When the new edge is connected to the generic triangle by two edges, these two edges and one edge of the generic triangle form a new triangle. Therefore, with probability p 2 , the generic triangle produces one child triangle. When the new edge is connected to the generic triangle by three edges, these edges and the edges of the generic triangle form three new triangles. Thus, with probability p 3 , the generic triangle produces three children triangles.
Any edge is called a type 2 object, and any triangle is called a type 3 object. We use subscript 2 for edges and subscript 3 for triangles. Thus, we denote by ξ i,j (t) the number of type j offspring of the type i generic object up to time t (i, j = 2, 3). Recall that ξ i,j , i, j = 2, 3, are point processes. Then gives the total number of offspring (that is both edges and triangles) of the generic edge up to time t. We can also see that is the number of all offspring (edges or triangles) of the generic triangle up to time t. We denote by τ 3 (1), τ 3 (2), . . . the birth times of the generic triangle, and we denote by ε 3 (1), ε 3 (2), . . . the corresponding total litter sizes. That is, at the ith birth event, the generic triangle bears ε 3 (i) children being either triangles or edges. The discrete random variables ε 3 (1), ε 3 (2), . . . are independent and identically distributed having distribution P(ε 3 (i) = j) = q j , j ≥ 1. By the above evolution process, we have We assume that the litter sizes are independent of the birth times. Let λ 3 be the life-length of the generic triangle. It is a finite, non-negative random variable. We assume that the reproduction terminates at the death of the individual. Therefore, ξ 3 (t) = ξ 3 (λ 3 ) for t > λ 3 . Then, the reproduction process of a triangle can be formulated as where Π 3 (t) is the Poisson process, S 3 (n) = ε 3 (1) + · · · + ε 3 (n) gives the total number of offspring of the generic triangle before the (n + 1)th birth event and by x ∧ y we denote the minimum of {x, y}. The survival function of the life-length. Let L 3 (t) denote the distribution function of the triangle's life-length λ 3 . Then, the survival function of λ 3 is where l 3 (t) is the hazard rate of the life-length λ 3 . We suppose that the hazard rate depends on the total number of offspring, so that with fixed positive constants b and c. Let λ 2 be the life-length of the generic edge. Then, ξ 2 (t) = ξ 2 (λ 2 ) for t > λ 2 . As the edge always gives birth to one offspring (which can be an edge or a triangle); therefore, is the total number of offspring of the generic edge, where Π 2 (t) is the Poisson process. We denote by L 2 (t) the distribution function of λ 2 . Then, the survival function of the life-length of an edge is where l 2 is the hazard rate of the life-length λ 2 . We suppose that l 2 is of the form l 2 (t) = b + cξ 2 (t).
We emphasize that we do not delete any edge or any triangle when it dies, because its ingredients can belong to other triangles or edges, too. Thus, dead triangles and edges will be considered as inactive objects not producing new offspring.
In Figure 1, an example is shown for our graph evolution model. For a clear view it contains only three birth steps after the initial time t = 0. The nodes of the ancestor are highlighted by red. The edges are labelled with the birth times t. The following objects appear in Figure 1

General Results
The survival functions.

Theorem 1. The survival function for a triangle is
The survival function for an edge is Proof. At the first part of the proof we omit subscripts 2 and 3, because the calculations are the same for edges and triangles. Let t > 0 and assume that Π(t) = k. Then, the first k birth events happened before time t. Thus, the birth times τ(1), τ(2), . . . , τ(k) and the corresponding litter sizes ε(1), ε(2), . . . , ε(k) are known. Therefore, the reproduction process ξ(u) is also known for u < t. By (5), a simple calculation shows that the survival function of an object is .
Let U * 1 , . . . , U * k be an ordered sample of size k from uniform distribution on [0, 1]. Then, the joint conditional distribution of the birth times τ(1), . . . , τ(k) given Π(t) = k, coincides with the distribution of tU * 1 , . . . , tU * k . Therefore because τ(i) = tU * i . The litter sizes ε(1), . . . , ε(k) are independent identically distributed random variables, which are independent also of U * 1 , . . . , U * k . Hence where we applied that U i is uniformly distributed. Using this and the total probability theorem, we find Therefore, the survival function for a triangle is Finally, the survival function for an edge is The mean offspring number. Let us denote by m i,j (t) = Eξ i,j (t) the expectation of the number of type j offspring of a type i mother until time t. where For any t ≥ 0, we have where where ε i,j (k) is the number of type j offspring of a type i mother at her kth birth event.

The Perron root and the Malthusian parameter.
Let be the matrix of the Laplace transforms. Direct calculation gives that the characteristic roots of M(κ) are The greater of the values 1 (κ) and 2 (κ) is called the Perron root, so is the Perron root. We assume that our process is supercritical; that is, For supercriticality, condition That value of κ for which the Perron root is equal to 1 is called the Malthusian parameter. Thus, using the usual notation in the theory of branching processes, α is the Malthusian parameter if (α) = 1. In this paper, we assume the existence of the Malthusian parameter. From relation (α) = 1 and (23), we obtain that the Malthusian α satisfies the equation Later, we use the eigenvectors of M(α). To this end, let α be the Malthusian parameter, and let (v 2 , v 3 ) be the right eigenvector of M(α) corresponding to eigenvalue 1 and Again, let α be the Malthusian parameter and let (u 2 , u 3 ) be the left eigenvector of M(α) satisfying condition u 2 v 2 + u 3 v 3 = 1. Direct calculation shows that

Asymptotic Theorems on the Number of Triangles and Edges
In this section, we use Proposition 4 from Section 8. So we should check the conditions given in Section 8. For condition (a) from Section 8, we should guarantee that not all measures m i,j are concentrated on a lattice. By Corollary 1, these measures are absolutely continuous, and thus it is satisfied.
Concerning condition (b1), we underline that we suppose the existence of a positive Malthusian parameter α. To this end, in this section, we assume that (25) has a finite positive solution α. We can check numerically the existence of this value. For (b2), we assume (24). Condition (c) from Section 8 will be checked later in the proofs of the results together with other conditions related to it. Now, we analyse condition (d). We can see from Corollary 1 that F(∞) and G(∞) are positive. Thus, we can concentrate on parameters r i and p i . If r 2 = p 1 = 0, then (d) is not satisfied; however, in this case, one can study separately the process of edges (it grows at any birth time by 1), and the process of triangles (this is described in [17]). If r 1 = 0 and p 2 + p 3 = 0, then (d) is not satisfied, and the evolution process is an alternating one. If either r 2 = 0 or p 1 = 0, then (d) is not satisfied.
To guarantee condition (d), in this section, we assume that 0 ≤ r 1 < 1, 0 < p 1 ≤ 1, and it is excluded that both r 1 = 0 and p 1 = 1 are satisfied at the same time. In this case, condition (d) from Section 8 is satisfied.
The denominator in the limit theorem. In the following theorem, we need the next formulae. In Section 8, we see that the denominator of m Φ ∞ in the limiting expression is independent of Φ, and it is It can be written in the form (and considering our two-dimensional case) Here, u i and v i are from Equations (26) and (27). Moreover, by Corollary 1 or by Proposition 1, we have that where Now, we turn to the number of edges and triangles. Recall that an edge is a type 2, and a triangle is a type 3 object. (24) is satisfied and (25) has a finite positive solution α. Assume that 0 ≤ r 1 < 1, 0 < p 1 ≤ 1 and it is excluded that both r 1 = 0 and p 1 = 1 are satisfied at the same time.

Theorem 2. Assume that
Let i E(t) denote the number of all edges being born up to time t if the ancestor of the population was a type i object, i = 2, 3. Then almost surely for i = 2, 3. Let iÊ (t) denote the number of all edges present at time t if the ancestor of the population was a type i object, i = 2, 3. Then almost surely for i = 2, 3. Let i T(t) denote the number of all triangles being born up to time t if the ancestor of the population was a type i object, i = 2, 3. Then almost surely for i = 2, 3. Let iT (t) denote the number of all triangles present at time t if the ancestor of the population was a type i object, i = 2, 3. Then, where and If i = 2, then ξ 2 (t) is the birth process of an edge, and the children can be both edges and triangles. Therefore, at each birth, there is one child. Therefore, where τ(1), τ(2), . . . are the jumps of the Poisson process Π 2 . In the Poisson process Π 2 (t) the distribution of the interarrival time (τ(i) − τ(i − 1)) is exponential with rate 1. Therefore, τ(i) has Γ-distribution Γ(i, 1). Using this, we have Let us denote by η i the interarrival time τ(i) − τ(i − 1). Let η 0 be an exponentially distributed random variable with rate 1 that is independent of M. Then, Therefore, the distribution of e −αη 0 (1 + M) coincides with the distribution of M. Therefore, using (40), we have From this, we find EM 2 = α + 2 2α 2 < ∞. Thus, (37) is true for i = 2. If i = 3, then ξ 3 (t) is the birth process of a triangle and the children can be both edges and triangles. Therefore, at each birth there are at most three children. Therefore, where τ(1), τ(2), . . . are the jumps of the Poisson process Π 3 . By the above calculation EM 2 < ∞, so (37) is true for i = 3.
If we show that ∞ 0 t 2 e −αt m i,j (dt) < ∞, for i, j = 2, 3, then conditions (c) and (iv) of Section 8 will be proved. Now, for i = 2 and j = 2, 3, we have from Corollary 1 To obtain (34), let Φ x (t) = 1 if x is an edge and it is present at t, and Φ To obtain (36), let Φ x (t) = 0 if x is an edge, and Φ x (t) = 1 if x is a triangle, and it is present at t. Therefore, EΦ 2 (t) = 0 and EΦ 3 Thus, (69) and (70) imply (36).

Generating Functions and the Probability of Extinction
The joint generating function of Π 2 (λ 2 ), ξ 22 (λ 2 ) and ξ 23 (λ 2 ). Recall that Π 2 is the Poisson process describing the reproduction times of the generic edge and λ 2 is its life length. Thus, is the joint distribution of the offspring size of the generic edge during its whole life and its last reproduction time. We have where τ i is the ith jumping time of the Poisson process Π 2 . Thus, it again shows that w i,j,k is the probability that the ith birth event is the last one that occurred before death, and the total numbers of the two types of offspring up to time τ i are equal to j and k, respectively. Now, consider the sequence Let ξ 2 (τ i−1 ) = m and assume for a while that τ i and τ i−1 are fixed. Then, using (4) and (5) for the hazard rate, we can calculate that, for fixed τ i and τ i−1 , We know that the increment (τ i − τ i−1 ) is exponential with parameter 1; therefore, At each birth step, the new individual can be either an edge or a triangle. Therefore, using the above calculations, the total probability theorem, and the independence of the type of the newly born individual and (Π 2 , λ 2 ), we have the following recursion for u i,j,k .
Now, by the definition of w i,j,k , we can see that where by (41),

b+c(j+k)
1+b+c(j+k) is the probability that the generic individual dies before the next birth event. Let . Then, from (42), we obtain the following recursion for the sequence v i,j,k where the initial values are v 0,0,0 = 1 1 + b and v 0,j,k = 0 for j = 0 or k = 0.
Now, we calculate the generating function G(x, y, z) of the sequence v i,j,k . We have First, multiplying with x i y j z k and then taking the sum of both sides of (43), we obtain where v 0,j,k , j = 0, 1, . . . is given by (44), and we define v i,j,k = 0 if j < 0 or k < 0. From this equation, we find Let h(t) = G(x, ty, tz). Now, substituting y with ty, z with tz in (45), we can obtain the following linear differential equation.
with the initial value condition Now, we can use the well-known method for linear differential equations. We obtain that the solution of the initial value problem (46) and (47) is We need the generating function of w i,j,k = v i,j,k (b + c(j + k)). It is From here, we obtain The joint generating function of Π 3 (λ 3 ), ξ 32 (λ 3 ) and ξ 33 (λ 3 ). Here, we study the offspring of a triangle. To distinguish the notation of this subsection and the previous subsection, but avoid too many subscripts, we use bar. Thus, here w i,j,k , u i,j,k , v i,j,k , G(x, y, z) and H(x, y, z) denote quantities relating offspring of the generic triangle. Recall that Π 3 is the Poisson process describing the reproduction times of the generic triangle and λ 3 is the life length of the triangle. Thus, is the joint distribution of the offspring size of the generic triangle during its whole life and its last reproduction time. We have where τ i is the ith jumping time of the Poisson process Π 3 . Thus, we again show that w i,j,k is the probability that the ith birth event is the last one that happened before death, and the total numbers of the two types of offspring up to time τ i are equal to j and k, respectively.
Let ξ 3 (τ i−1 ) = m, and assume for a while that τ i and τ i−1 are fixed. Then, using (4) and (5) for the hazard rate, we can calculate that, for fixed τ i and τ i−1 , We know that the increment (τ i − τ i−1 ) is exponential with parameter 1; therefore, At each birth step, the new individual can be either an edge or a triangle. Therefore, using the above calculations, the total probability theorem, and the independence of the type of the newly born individual and (Π 3 , λ 3 ), we have the following recursion for u i,j,k .
. (52) Now, by the definition of w i,j,k , we can see that where by (51),

b+c(j+k)
1+b+c(j+k) is the probability that the generic individual dies before the next birth event.
. Then, from (52), we obtain the following recursion for the sequence v i,j,k where the initial values are v 0,0,0 = 1 1 + b and v 0,j,k = 0 for j = 0 or k = 0.
Now, we calculate the generating function G(x, y, z) of the sequence v i,j,k . We have First, multiplying with x i y j z k and then taking the sum of both sides of (53), we obtain where v 0,j,k , j = 0, 1, . . . is given by (54) and we define v i,j,k = 0 if j < 0 or k < 0. From this equation, we find (1 + b) G(x, y, z) − 1 1 + b + ycG y (x, y, z) + zcG z (x, y, z) = = p 1 xyG(x, y, z) + p 2 xzG(x, y, z) + p 3 xz 3 G(x, y, z). (55) Let h(t) = G(x, ty, tz). Now, substituting y with ty, z with tz in (55), we can obtain the following linear differential equation.
with the initial value condition One can see that the solution of the initial value problem (56) and (57) is With t = 1, we obtain that Therefore, the generating function of From here, we obtain Actually, m i,j (∞) is the expected offspring number of the embedded Galton-Watson process.
Let s 2 and s 3 denote the probability of extinction of our process when the ancestor is an edge, resp. triangle. Theorem 3. Assume that 0 ≤ r 1 < 1, 0 < p 1 ≤ 1 and it is excluded that both r 1 = 0 and p 1 = 1 are satisfied at the same time. Let be the Perron-Frobenius root of M. If ≤ 1, then s 2 = s 3 = 1. If > 1, then s 2 < 1 and s 3 < 1. In any case, (s 2 , s 3 ) is the smallest non-negative solution of the vector equation (s 2 , where f 2 and f 3 are given in Corollaries 2 and 3. Proof. We apply Theorem 7.1 in Chapter 1 of [19]. By Corollary 1, m i,j (0) = 0 and m i,j (t) is finite for any i, j. Therefore, by Theorem 7.1 in Chapter 3 of [19], the extinction of our original process has the same probability as the extinction of the embedded Galton-Watson process. Thus, we can apply Theorem 7.1 in Chapter 1 of [19]. Here, M is the matrix of the expected offspring numbers of the embedded Galton-Watson process. Now, M is positively regular because we assume that 0 ≤ r 1 < 1, 0 < p 1 ≤ 1 and it is excluded, that both r 1 = 0 and p 1 = 1 are satisfied at the same time. Thus, our result follows from Theorem 7.1 in Chapter 1 of [19].

The Asymptotic Behaviour of the Degree of a Fixed Vertex
The process of the 'good children'. To describe the degree of a fixed vertex, we introduce a new branching process that we call the process of 'good children'. This process contains those objects that contribute to the degree of the fixed vertex. We can see that a newly born vertex can have 1 or 2 edges if its parent is an edge object and 1, 2 or 3 edges if its parent is a triangle object.
First, we consider the case when the newly born vertex has one edge, and thus, at the beginning, it belongs to an edge object. In this paragraph, we call this edge the 'parent' edge. We fix the newly born vertex. Then, we distinguish those children objects of the 'parent' edge, which contribute to the degree of our fixed vertex. We call a child object of the 'parent' edge a 'good child' if it contains our fixed vertex. We can see that only the 'good children' and their 'good children' offspring can contribute to the degree of the fixed vertex. Then, the distribution of the number of 'good children' at a reproduction event of the 'parent' edge is whereε 22 denotes the number of edge type 'good children' andε 23 denotes the triangle type 'good children'. We have to consider the reproduction process of the 'good child', which is the following Limit results for the degree. We have already mentioned that the 'good children' and only they can contribute to the degree of the fixed vertex. Thus, its degree is equal to the initial degree plus the number of 'good children'. Let 2C (t) be the degree of a fixed vertex at time t after its birth in the case when the vertex belongs to an edge at its birth. Similarly, 3C (t) is its degree in the case when the vertex belongs to triangle at its birth. Up to an additive constant, iC (t) is the number of 'good children' offspring of an i type 'parent' object at time t. It is the sum of the number of edge type 'good children' iẼ (t) and the triangle type 'good children' iT (t). To apply Proposition 4, we can use the same method as in Theorem 2. Thus, for the edges, we can again use the random characteristic Φ x (t) = 1 if x is an edge and Φ x (t) = 0 if x is a triangle, but the underlying process is the process of 'good children'. This is similar for triangles. Therefore, we have almost surely for i = 2, 3, where 2W and 3W are positive on the event of non-extinction of the 'good children'. The last case is when the newly born vertex has three edges. Then, three triangles contribute to the degree of that vertex. Let 3C (t) be the degree of this vertex. Then, 3C (t) is the sum if 'good' offspring of three triangles. Thus, almost surely, where 3W1 , 3W2 , 3W3 are independent copies of 3W . Checking the conditions of Proposition 4 for the 'good children' process. To complete the previous reasoning, we should check the conditions of Proposition 4. First, we find the the denominator in the limit theorem that is we calculateD. By Section 8, we see thatD Here,ũ i andṽ i are the eigenvectors. Moreover, whereα is the Malthusian parameter in the process of 'good children' and A , B denotes the derivatives given in (31) and (32). Condition (a) of Proposition 4 is true because the measuresm i,j are non-lattice as they are absolutely continuous. For condition (b1), we assume the existence of a positive Malthusian parameter. That is, we assume that (64) has a finite and positive solutionα. Condition (b2) is true, because we assume that˜ (0) > 1. Condition (c) is a consequence of Section 4, becausem i,j (t) has shape cm i,j , where c is positive number.
To guarantee condition (d), in this section, we assume that 0 ≤ r 1 < 1, 0 < p 1 ≤ 1, and it is excluded that both r 1 = 0 and p 1 = 1 are satisfied at the same time. Conditions (i)-(ii)-(iii) and (v) are true because of the shape of Φ. Conditions (iv) and (vi) are consequences of ξ i,j (t) ≤ ξ i,j (t) as one can see from the proof of Theorem 2.  For each parameter set, we stored the mentioned measurements only in integer time steps, and then we took the average of 100 simulated processes. In Figure 4, an example is shown for a specific parameter set (r 1 = 0.1, p 1 = 0.2, p 3 = 0.6, b = 0.25, c = 0.25). The values of the averages are plotted by dots. In each case, we fitted a regression line (plotted by continuous red line) to the last 9 values. We can see that the fit is perfect, thus, supporting our theorem.
Our main goal was to obtain a 95% confidence interval for the slope of the linear regression line, as that was our simulated approximation of the Malthusian parameter α. Table 1 contains the boundaries of the 95% confidence intervals for α. The columns labelled with 2.5% and 97.5% refer to the lower and the upper bounds obtained from simulations, while the column ofα refers to the numerical solution of Equation (25).
For each fixed parameter set {r 1 , p 1 , p 2 , b, c}, we present the confidence intervals calculated from the number of edges being born (E) resp. being present (Ẽ) and from the number of triangles being born (T) resp. being present (T) up to time t = 14. The confidence intervals containing the numerical Malthusian parameterα are highlighted with the * symbol. We see that any confidence interval is narrow, and it either containŝ α, orα is very close to the interval. These results show that the approximation is good for moderate values of t.  Finally, we present some simulation results for Theorem 3, that is, for the probability of extinction of the evolution process. We made the following computer experiment for any fixed parameter set {r 1 , p 1 , p 2 , b, c} and for type 2 and type 3 ancestors. We started to generate the process. If this process reached 2 10 birth steps, then we stopped it and considered it as a non-extinct process. Otherwise, when the process did not reach 2 10 birth steps, then the process died out. Applying the above method, we generated 10 5 processes for each parameter sets and counted the relative frequencies of the processes being extinct.
In Table 2, we show some of the results. Column Ancestor contains the type of the ancestor. In the column Numeric we show the numeric solution of the non-linear equation in Theorem 3. We used Julia's trust region method. Column Simulation contains the relative frequencies extracted from the simulations. The simulation results slightly underestimate the numeric values. This is reasonable because we stopped all processes at a fixed time.

Basic Facts on Branching Processes
In our paper, we use known results of the theory of continuous-time branching processes. The single type general Crump-Mode-Jagers branching processes have been described e.g., in [21][22][23]. The general multi-type branching processes have been studied, e.g., in [11,19,24].
Here, we give a short description of the general multi-type branching processes based on [11]. The individuals of this process can be of p different types, which we denote by 1, 2, . . . , p. Any individual x is described by the quantities λ x , ξ x , Φ x , Ψ x , . . . . The quantities λ x , ξ x , Φ x , Ψ x , . . . are independent copies of the quantities λ, ξ, Φ, Ψ, . . . . Thus, we should give the definition of λ, ξ, Φ, Ψ, . . . , which we consider as the quantities corresponding to the generic individual.
The lifetime λ is a non-negative random variable which is not necessarily independent from the reproduction. The lifetime distribution is L(t) = P(λ ≤ t). The reproduction process is ξ i (t) = ξ i,1 (t), . . . , ξ i,p (t) , t ≥ 0. Here, the random point process ξ i,j describes the births of type j offspring of a type i mother. ξ i,j (t) gives the number of type j offspring of a type i mother up to time t. ξ i,j is determined by the birth events and the numbers of offspring. The process starts at time t = 0 with one individual called the ancestor and denoted by x 0 . When a child is born, it starts its own reproduction process and so on. The birth time of the individual x is denoted by σ x .
Let Φ(t) be a non-negative random function that describes a certain aspect of the life history of the individual. It is usually assumed that Φ(t) = 0 for t ≤ 0. Then, Φ(t) is called a random characteristic. Let Ψ(t) be another random characteristic. Thus, the behaviour of the individual x is described by ξ x , λ x , Φ x , Ψ x , . . . .

for any ancestor.
Then, lim t→∞ e −αt is likely, where i is the type of x 0 , x 0 Y ∞ is an a.s. non-negative random variable depending on the type of the ancestor x 0 but not depending on the choice of Φ.
If, in addition, we assume that then E( x 0 Y ∞ ) = 1, x 0 Y ∞ is positive with positive probability, and x 0 Y ∞ is a.s. positive on the survival set.
The proof is a simple consequence of Theorem 2.4 and Proposition 4.1 of [11].

Discussion
In this paper, a new network evolution model was introduced. This model was inspired by those networks where small substructures play important role. In social life, such substructures could be a group of friends. In the theory of networks, these substructures are called motifs. In this paper, for the sake of simplicity, we consider only two types of substructures, the edges and the triangles. The novelty of the paper is the usage of a two-type continuous time branching process to describe these two types of interactions. Thus, despite [7,10], the theory of multi-type branching processes was applied for certain substructures of the network and not just for the nodes. Our paper extends the former studies of [16,17], where only one type of interaction was considered.
In this paper, we proved that the magnitude of the number of triangles on the event of non-extinction is e αt , where α is the Malthusian parameter. We obtained similar results for the number of edges. We also studied the degree process of a fixed vertex and the probability of the extinction. Our results are similar to the ones obtained for the simpler models in [16,17]. In addition to mathematical proofs, the results were illustrated by simulations.
In future extensions of the model, more than two types of substructures can be studied using the theory of multi-type branching processes.