Construction of k-Variate Survival Functions with Emphasis on the Case k = 3

The purpose of this paper is to present a general universal formula for k-variate survival functions for arbitrary k = 2, 3, ..., given all the univariate marginal survival 
functions. This universal form of k-variate probability distributions was obtained by means of “dependence 
functions” named “joiners” in the text. These joiners determine all the 
involved stochastic dependencies between the underlying random variables. However, 
in order that the presented formula (the form) represents a legitimate survival 
function, some necessary and sufficient conditions for the joiners had to be 
found. Basically, finding those conditions is the main task of this paper. This 
task was successfully performed for the case k = 2 and the main results for the case k = 3 were formulated as Theorem 1 and Theorem 2 in Section 4. Nevertheless, 
the hypothetical conditions valid for the general k ≥ 4 case were also formulated in Section 3 as the (very convincing) 
Hypothesis. As for the sufficient conditions for both the k = 3 and k ≥ 4 cases, the full generality was not achieved since two restrictions 
were imposed. Firstly, we limited ourselves to the, defined in the text, 
“continuous cases” (when the corresponding joint density exists and is 
continuous), and secondly we consider positive stochastic dependencies only. 
Nevertheless, the class of the k-variate distributions which can be constructed is very wide. The 
presented method of construction by means of joiners can be considered 
competitive to the copula methodology. As it is suggested in the paper the 
possibility of building a common theory of both copulae and joiners is quite 
possible, and the joiners may play the role of tools within the theory of 
copulae, and vice versa copulae may, for example, be used for finding proper 
joiners. Another independent feature of the joiners methodology is the 
possibility of constructing many new stochastic processes including stationary 
and Markovian.


Introduction
This work includes a continuation of our previous papers [1] [2] [3] and [4] on analysis and construction methods of multivariate survival functions, given their univariate marginals.
This subject was widely developed in the literature since the late 1950s and 60s (for example, see [5] [6] [7] [8] [9]) until nowadays. Many methods of construction of bivariate and multivariate probability distributions were developed more recently, only to mention the conditionally determined bivariate and multivariate distributions in [10] [11] [12]. Some similar but different methods of constructing models, under the name "parameter dependence method", can be found in [13] [14]. As it turned out many of the multivariate distributions obtained by that method can also be obtained by the "method of triangular transformations" (especially by pseudoaffine and pseudopower transformations), see for example [15].
Many other options for construction one can be found in the references of [7]. The approach that is developed here and in [1] [2] [3] has its origin [4] in the model which is the Aalen version [16] of the famous Cox model [17] for stochastic dependence.
However, as a result of further development we finally were able to formulate (here and in [1] [2] [3]) a version of an emerging theory independent of its Aalen origin.
Thus, in our approach, as presented here and in [1] [2] [3], both the construction and the universal characterization of multivariate models rely on incorporating the so called "joiners". These joiners we propose to name "Aalen factors" (but we will use the name "joiner" for short). They are the functions which fully determine all underlying stochastic dependencies between the considered random variables.
As a result, the so obtained k-variate survival functions gain a nice factored form ( ) ( ) ( ) ( ) ( ) , , ,  is a valid k-variate survival function.
This task formulation, in terms of the survival functions, resembles the task of finding a proper copula for given univariate cdfs, say, in order to obtain a k-variate cdf. as a stochastic model for some investigated reality (i.e., "outside of mathematics") represented by a given set of data.
It becomes then clear that, given, "essentially the same" input data, , the "method of joiners" stands as an alternative (and competitive [18]) method to the copula methodology.
On the other hand, as we point out in Section 2, there is a possibility to formulate a common general theory of both copula and joiner representations of bivariate and, possibly, also k-variate (k ≥ 3) probability distributions.
This possibility follows from the fact that both representations (by copulas and by joiners) are equivalent as describing the same probability distribution.
Therefore, any copula uniquely determines a corresponding joiner, and any joiner uniquely determines a copula. This fact indicates the way to find more copulas through joiners, as well as more joiners that correspond to known copulas. Moreover, the method of joiners may, likely, be used to develop more copulas theory for higher dimensions.
These remarks only signalize the possibility of such a common theory. We intend to develop this in more detail in the nearest future.
At the moment, we rather concentrate on the joiner representation and the joiner based methods of models construction.
In the next section we provide a short introduction to the joiners' theory by providing a slightly different (as compared with our previous works) formulation of bivariate case. This case is fundamental to our considerations in sections 3 and 4 since for the k-variate (k ≥ 3) survival functions we restrict our attention to "bi-dependence" only, which means only bivariate joiners may be different than 1 (see, [2]). This restriction dramatically simplifies theory of k-variate distributions as compared with the general theory (for arbitrary k) developed in [2]. In Section 3 we provide only a general scheme of k-variate distributions' analytical form (for arbitrary k) under the bi-dependence assumption.
As for the k-variate distribution for which the joiner (in this case the joineris reduced to aproduct of bivariate joiners that correspond to all, or to some only, bivariate marginals) must fit the given k univariate marginals, we only formulate the main result as Hypothesis.
Even though the Hypothesis is very convincing, we were not able to provide a formal proof. The arguments, which made us strongly believe it holds, mostly (but not exclusively) follow from the fact that the same pattern as we presume holds for an arbitrary k, holds for the cases k = 2 and k = 3.
As mentioned, the (hypothetical) conclusion, extending these two cases to all k, is not only based on that analogy, but also on the naturalness of the assertion.
The case k = 3 is the subject of Section 4. The main thesis about the 3-dimensional model is formulated and proven there in two theorems.
The general case, i.e., k-variate distribution for any 1, 2, k =  , actually defines an infinite sequence of random variables 1 2 , , , , k X X X   which, as it is pointed out at end of Section 3, satisfies the Kolmogorov consistency condition.
So, we actually defined a class of discrete time (now, k represents "time") stochastic processes with a variety of interesting special cases.
Clearly, such processes need not to be very complicated if we assume that most of the bivariate joiners reduce to 1. Based on such possibility we may construct m-Markovian processes for 1, 2, m =  (they are Markovian if m = 1).
Also by adopting natural assumptions we may gain stationarity of some constructed stochastic processes. It is an exciting possibility that having only one bivariate distribution of any neighboring random variables, say 1 , k k X X − , we may easily construct both stationary and Markovian stochastic processes.
This subject is only touched upon in Section 3, but is out of scope of this work. At the end of this Introduction we must notice that, according to our common assumption that every joiner is not bigger than 1 (i.e., its corresponding representation, by the defined throughout the text function , is nonnegative) we limited ourselves, to positive stochastic dependencies. Extension to models also comprising negative dependencies is quite possible, but requires a more complex theory.

The Bivariate Case
Before constructing new classes of multivariate and, especially, tri-variate survival functions, in this section we repeat and slightly modify the bivariate cases which involve their bivariate universal forms, referring for more details to our previous papers [1] [3] [4].
Thus, according to those papers, any bivariate survival function ( ) ( ) of an arbitrary random vector ( ) 1 2 , X X can be represented as: is a continuous function satisfying (4), equality (3) may be rewritten into the form The question whether, for the marginals given by the hazard rates ( ) x λ , the function ( ) x x ∂ ∂ ∂ of ( ) 1 2 , R x x , as given by (5), are equal to each other and are nonnegative for all 1 x , 2 x .
Obviously they do exist and are continues by the earlier assumed continuity of the functions ( ) The nonnegativity requirement for ( ) is equivalent to the simple common fact that the joint density of any random vector ( ) 1 2 , X X , if it exists, must be nonnegative.
As it follows from the form of (3), other properties, characterizing that density's antiderivative (the cdf.) are satisfied automatically.
Thus, to obtain the sufficient condition for ( ) ( ) where ( ) 1 2 , R x x is given by (5).
After calculating the second mixed derivative from (5), then simplifying underlying expressions and setting expressions with negative signs to the right-hand side of the inequality, one obtains (6) in the form of the following integral inequality: Now, the task of finding the bivariate distribution, given the marginals as represented by the hazard rates x λ , reduces to solving integral inequality (7) with respect to the only unknown function . This means, any solution x λ . So, in the continuous case, the set of all the solutions ( ) (7) uniquely determines the set of all bivariate distributions having the same fixed marginals. All the functions ( ) , satisfying conditions (4) and (7), will be denoted by the symbol . So that whenever writing ( ) we will mean the function representing a proper joiner and the corresponding, by (3), function denoted by " ( ) 1 2 , S x x " will be treated as a legitimate bivariate survival function.
In the case ( ) if the inequality: holds then (7) holds too. Thus, the conditions ( ) Notice, however, that condition (8) is not a necessary condition. Nevertheless, any solution of (8) is a solution of (7) and is very easy to find. The simplest set of such solutions obviously is given by , where for the constant parameter a (to be statistically estimated) we require 0 1 a ≤ ≤ .
One then obtains as a special case of (3) the following model: Model (9), which is somehow related to the first Gumbel bivariate exponential [6], is the most natural (and, possibly, kind of "canonical") solution to the problem of finding the joint distribution of random variables X 1 , X 2 , given the margins represented by the hazard rates.
We expect many applications of this model according to the common opinion that the simplest ("but not yet simpler") models mostly are the best reflections of modeled realities.
Model (9) can be generalized to the following one: where an additional factor of the middle term in the exponent of(10) is any continuous function ( ) 1 2 , c x x satisfying: In particular, we propose the following natural model, with where the constant real parameter γ (to be estimated) satisfies 0 γ ≥ . [18] to the theory of copulas [9]. However, it can be shown that the two theories are equivalent in the sense that there is a one-to-one relationship between joiners and copulas, at least in the bivariate case. Thus, having a joiner one immediately obtains the corresponding unique copula and vice versa. As it is well-known, every copula is "good" to any pair of marginal distributions, but, as follows from (7), not every "joiner" fits the marginals given by hazard rates

The theory of joiner representations (as developed in this work and in [1] [2] [3] [4]) is competitive
. On the other hand substituting into any copula the given marginals, one obtains back (through that copula) the joiner that fits the substituted marginals. Thus, in such a way, one can obtain all the proper joiners through all the copulas that are known.
It is important to notice that, in applications to practical problems, it is easier to find a joiner that reflects a modeled reality than a proper copula (That "easiness" in finding a proper joiner follows from the fact that the models involving joiners are closely related to the Aalen [16] version of the Cox [17] model for stochastic dependencies.).
On the other hand, such a proper joiner determines the corresponding copula. This facilitates the choice of proper copula. In that sense the joiner approach may be considered as an enrichment of the copula methodology. This subject of a possible common theory of copulas and joiners will be included in our next publication which is now in preparation.

A Class of k-Variate Survival Functions
Before starting a more detailed analysis of tri-variate survival functions, which is our main goal, we first introduce a simplified version of k-variate survival functions for any k ≥ 3. This version turns out to be a special case of the most general k-variate model presented by formula (1) in [2].
The model here presented is much simpler as it only describes the "pairwise stochastic dependence" which was defined in [2]. According to the terminology in [2] such distributions are "three-independent".
As it will be seen, this simplified case still preserves a significant amount of generality.
Consider the following simplified formula for a k-variate survival function of the random vector ( ) As mentioned, this formula may be considered as a special case of formula (1) in [2]. Now, however, as we assume pairwise dependence only, the joiner in (12) reduces to the following product of bivariate joiners: Comparing (13) to the joiner as present in formula (3) of [2], one sees that (13) is obtainable from the most general case considered in [2] by setting to 1 all the joiners which are not bivariate.
Resuming, formula (12) can be rewritten as: where any single bivariate joiner In the continuous case this means that all the k hazard rates Naturally, for all pairs of hazard rates ( ) ( ) hold for all , i j x x , given any fixed pair of indexes i, j.
Inequalities (17) have the same structure as inequality (7). With the above continuity assumptions, formula (14) can be rewritten in the following exponential form: Notice that if we again substitute into (18) which upon condition (17) is a well-defined bivariate survival function. Remark It's easy to find out that if we set in (18) any k r − (1 r k ≤ < ) variables (among all the variables 1 , , k x x  ) to 0, we obtain the r-dimensional marginal survival function of the remaining random variables, say 1 , , i ir X X  which, syntactically, has an identical form as the k-dimensional survival function (18). This observation implies the following: 1) Since the pattern given by (18) is valid for every 2,3, k =  we, in fact, defined (at least theoretically) an infinite sequence of probability distributions.
2) The above observation on r-dimensional marginals of every k-dimensional distribution, (for each 1, 2, , 1 r k = −  ) clearly indicates that for the underlying sequence of random variables 1 , , , k X X   ( 1, 2, k =  ) the Kolmogorov and Daniels consistency requirement is satisfied in this case. Therefore a fairly wide class of stochastic processes { } 1,2, k k X =  with discrete "time" k is defined.
3) These stochastic processes will be significantly simplified to Markovian if we set all the joiners Thus, under the foregoing assumption, only the "adjacent" random variables, say  . Of course, in applications, the choice of joiner (i.e., which one was to be set to 1 ) depends on the character of a modeled reality. Anyway, there is no need always to rely on Markovian (m = 1) processes. 4) Suppose, for the defined above Markovian (m = 1) stochastic processes, we assume that all the underlying hazard rates have the same functional form, i.e., Also assume one functional form for all the bivariate joiners present in defining formula (14).
Then the so defined Markovian (as well as any m-Markovian) process will be stationary, and all we need to analyze and for any further computations regarding the whole so encountered process, it is enough to know one bivariate survival function, say for all 2,3, i =  . In more general cases (not necessarily "continuous") one may consider, instead, the survival function in the form: 5) As one can realize, an interesting new theory of both random vectors and stochastic processes emerges. This is, however, not in the scope of this work and we postpone that subject for future research. □  . Nevertheless, we have reasons to propose the following hypothesis which, unfortunately, at the moment, we are not able to prove rigorously.

As for (19) we already know that if [for all pairs ( )
In the general k-variate case an eventual rigorous proof would, possibly, re-quire to establish some common pattern that would comprise the k th mixed de- from the right-hand side of (18). Such a pattern should be valid for any 3, 4, k =  , and it is more than finding the k th derivative for any particular relatively small k.
As for the hypothesis we presume what follows: Hypothesis: For any given univariate marginals of a random vector ( ) represented by the hazard rates consider the following expression: where for any pair of indexes i, j such that 1 i j k ≤ < ≤ the continuous func- Moreover, let the nonnegative real constants ij c satisfy: Then expression (21) represents a valid k-dimensional survival function of ( ) . □ If the above hypothesis holds, then from (22) one obtains the following specific k-variate model: where the nonnegative constants ij c satisfy (23).
Model (24) is an extension of the similar bivariate model (9) as well as the tri-variate model (26) given in the next section.
If the Hypothesis holds, then we also can generalize (24) to the following: for any set of continuous functions For example, one may choose  with a given set of real nonnegative constants ij b     . However, if the hypothesis holds, model (24) seems to be the most natural in the class of k-dimensional survival functions in the continuous case. We expect many applications of (24) in multivariate survival analysis.
Although we do not have at our disposal a formal proof of the considered Hypothesis, an argument that may support it is that it holds, in particular, for k = 2 (Section 2) and for k = 3 (next section). So the truthfulness of the hypothesis (for all k) is (unfortunately) based solely on that analogy. Nonetheless, we are strongly convinced it holds in all cases k ≥ 2. In a case of occurrence an essential difficulty in finding a formal proof for dimensions higher than 3, in applications (only), some statistical arguments for the cases 4,5, k =  may possibly be applied. Another way out might be the use of CAS (Computer Algebra Systems) such as MAPLE or MATEMATICA for underlying computations.
In the case when possessing a k-variate "model" (k ≥ 4) such as (21) is especially important for a given practical purpose one might, eventually, take a (slight) risk and adopt some model (21) together with condition (23), and then try to verify it statistically. This approach may turn out to be useful, at least from a practical point of view. However, such "purely statistical" arguments would not be equivalent to a proper analytical (mathematical) proof.

The Continuous "Canonical" Case
Consider the following 3-dimensional version of (24): where we assume that all three continuous hazard rates ( ) ( ) ( ) are never zero. This assumption may, possibly, be weakened. We adopt it only for simplification of our calculations.
The question now to be answered is, for which hazard rates and for which values of the constant parameters 12 13 23 , , c c c does expression (26) represent a valid 3-variate survival function.
The answer to this question, together with proper restrictions, we formulate as the following: Theorem 1. Given is a random vector ( ) , , X X X whose univariate marginal survival functions are represented by any given continuous hazard rates , which never are zero.
The function ( ) , , X X X , if the nonnegative coefficients ij c (1 3 i j ≤ < ≤ ) in (26) satisfy: 12 are satisfied by elementary properties of Riemann integrals.
Therefore, we may return to our original notation for the joiners, i.e., . Now consider the sufficient condition for (26) to be a survival function. Since for the, here considered, "continuous case" the derivative ( ) exists and is continuous, the sufficient condition reduces to the following inequality ( ) ( ) which must hold for all the triples ( ) , , where all the expressions A, B, C and all the constants 12 13 23 , , c c c are nonnegative. Also, A, B, C are increasing functions of the variables 1 2 3 , , x x x , while they all are 0 at ( ) ( ) 1 2 3 , , 0, 0, 0 x x x = . It is easy to solve inequality (30) since it is equivalent to inequality (31). (31) is quite easy to solve and get to the conclusion that when it holds, (30) does too whenever the nonnegative constants ij c satisfy: 12 13 23 The very important conclusion from the above form (31) of (30) and its set of solutions (given by (32)) is that both inequalities' truthfulness does not depend on either values or functional forms of the involved three continuous hazard rates ( ) ( ) ( )  At the end of this subsection notice that setting any of the variables 1 x , 2 x , For example setting 3 0 x = in (26) one obtains (9), with 12 c a = .

A More General Continuous Case
One of the most general forms of "continuous" 3-dimensional survival functions, which also comprise the case considered in Section 4.1, can be defined as follows: , which, as in the case of the models (26), restricts the considerations to positive [but, possibly, to all positive] stochastic dependencies only. As a matter of fact that nonnegativity assumption is not necessary, and is adopted only for simplicity reasons.
We will prove the following theorem: Theorem 2 Given that (32) and (34) hold, formula (33) defines a class of valid 3-variate survival functions.
Proof. The argumentation is mainly based on the already proven validity of (26) as defining survival functions. As in the case (26), we need to check if the following inequality (similar to (29)) holds: After similar computations as in the previous case we arrive at the following inequality equivalent to (35): where from the first inequalities in (34) [for all 1 3 i j ≤ < ≤ ] follows the nonnegativity of * * * , , A B C . Now, upon dividing both sides of (37) by the (always positive) product ( ) ( ) ( ) The equivalence of (40) with (36) and (37) does not depend on the always nonnegative values of the, given by (39), expressions for * * * , , A B C , even if they stricted to investigation of distributions of k-dimensional random vectors only.
As it was pointed out in that Remark, the method can as well be applied to the construction of stochastic processes with discrete time (the Kolmogorov-Daniels consistency conditions are always satisfied). The so constructed processes may possess very nice properties such as m-Markovianity (the Markovianity when m = 1) and stationarity. Easiness of such constructions indicates the significant value of the created "joiners theory".