Generalization of the Change of Variables Formula

with Applications to Residual Flows

Niklas KoenenuniBremen,bips
Marvin N. WrightuniBremen,bips
Peter MaaßuniBremen
Jens BehrmannuniBremen

uniBremenFaculty of Mathematics and Computer Science, University of Bremen, Bremen, Germany
bipsLeibniz Institute for Prevention Research and Epidemiology – BIPS, Bremen, Germany

Niklas<EMAIL_ADDRESS>
Machine Learning, ICML, Change of Variables, Generalization of Flows, Residual Flows, Measure Theory

Normalizing flows leverage the Change of Variables Formula (CVF) to define flexible density models. Yet, the requirement of smooth transformations (diffeomorphisms) in the CVF poses a significant challenge in the construction of these models. To enlarge the design space of flows, we introduce $\L$-diffeomorphisms as generalized transformations which may violate these requirements on zero Lebesgue-measure sets. This relaxation allows e.g. the use of non-smooth activation functions such as ReLU. Finally, we apply the obtained results to planar, radial, and contractive residual flows.

§ INTRODUCTION

The term normalizing flow refers to a concatenation of arbitrarily many simple transformations such that together they describe a transformation of desired flexibility and expressiveness. Formally, a transformation $f:Z \to X$ denotes a diffeomorphism, i.e., a bijective mapping where both $f$ and $f^{\minus1}$ are continuously differentiable. The crucial reason why transformations are considered in normalizing flows is the validity of the Change of Variables Formula (CVF) for a probability density $p_Z$ on $Z$ described by

\begin{align}\label{Original_CVF}
\int_{f^{\minus1}(A)} p_Z(z)\ d\lambda^n(z) = \int_{A} p_f(x)\ d\lambda^n(x),
\end{align}

where $p_f(x):= p_Z(f^{\minus1}(x)) |\det J_{f^{\minus1}}(x)|$ and $A\subseteq X$. This formula provides an explicit expression of the density $p_f$ induced by $f$ on the target space $X$, which includes the determinant of the Jacobian as a volume correction term.

Based on this definition of a transformation, however, it is generally not accurate to use non-smooth activations, such as ReLU, Leaky ReLU, or ELU with $\alpha \neq 1$, in the design of normalizing flows. These usually cause the flow to become non-differentiable on a set with no volume w.r.t. the Lebesgue measure $\lambda^n$, hence no diffeomorphism. In measure theory, these sets are called $\lambda^n$-null sets or only null sets for short and are negligible in integration.

We demonstrate that the requirements for a flow can be significantly weakened by excluding null sets from the base and target space while preserving the validity of the CVF. There are remarks on using almost everywhere (a.e.) differentiable activation functions in <cit.> or <cit.>, yet both works lack a proof of the validity of such transformations to define flows. In our work, we provide such proofs for an even more general statement. At the same time, we discuss the probabilistic background of normalizing flows to induce a well-defined density in the end. Furthermore, we point out that not every generalization of the CVF found in the mathematical literature is immediately suitable for flows. Finally, we put a special emphasis on the applications to residual flows. In doing so, we prove that non-smooth activations are also valid for both planar and radial flows <cit.>, as well as for contractive residual flows <cit.>.

§ BACKGROUND ON CVF IN PROBABILITY THEORY

The basic idea behind normalizing flows is to transform a known and tractable probability space into a more complex one. Mathematically, a probability space $(Z, \mathcal{A}_Z, \mathbb{P}_Z)$ is composed of a set $Z$ equipped with a $\sigma$-algebra and a probability measure $\mathbb{P}_Z$. In the target space $(X, \mathcal{A}_X, \mathbb{P}_\text{data})$, only the set and $\sigma$-algebra are fixed, and the data distribution $\mathbb{P}_\text{data}$ is unknown. For simplicity, we only consider open subsets of $\R^n$ and trace $\sigma$-algebras of the Lebesgue algebra $\L$ in the following, i.e., $\mathcal{A}_Z = \L(Z)$ and $\mathcal{A}_X = \L(X)$. The trace $\sigma$-Algebra is a restricted $\sigma$-Algebra on a subset defined by $\L(Z):= \{A \cap Z \mid A \in \L\}$. Besides, we assume that the distribution $\mathbb{P}_Z$ is absolutely continuous w.r.t. the $n$-dimensional Lebesgue measure $\lambda^n$, i.e., $\lambda^n$-null sets have a $\mathbb{P}_Z$-probability of zero. Therefore, the existence of a probability density $p_Z$ follows by Radon-Nikodym's theorem (). For more information on measure theory, see <cit.> or <cit.>.

At the lowest level, a transformation $f:Z \to X$ has to be at least an $\mathcal{A}_Z$-$\mathcal{A}_X$-measurable mapping (i.e., a random variable) in order to induce a distribution on the target space by the so-called pushforward measure:

\begin{align*}
\mathbb{P}_f(A) := \mathbb{P}_Z\left(f^{\minus1}(A)\right) && (A \in \mathcal{A}_X). %\label{Background:PushforwardMeasure}
\end{align*}

Under the assumption that the base distribution $\mathbb{P}_Z$ has a density function $p_Z$, there also exists an integral representation for the pushforward measure, i.e.,

\begin{align}
\mathbb{P}_f(A) = \int_{f^{\minus1}(A)} p_Z(z)\ d\lambda^n(z). \label{Background:PushforwardMeasureWithDensity}
\end{align}

[Transformations for the CVF]
The function $f:Z \to X$ between two open sets $Z,X \subseteq \R^n$ is a diffeomorphism; or equivalently expressed by the inverse function theorem, $f$ is bijective, continuously differentiable and without critical points.

We note that $z \in Z$ is a critical point if the Jacobian-determinant vanishes in this point, i.e., $\det J_f(z) = 0$. In particular, a critical point $z$ indicates that the inverse is non-differentiable or non-continuously differentiable in $f(z)$.

If the mapping $f$ satisfies Assumption <ref>, the CVF from eq. (<ref>) holds, and we can extend the expression (<ref>) of the distribution $\mathbb{P}_f$ by

\begin{align}
\mathbb{P}_f(A) = \int_{f^{\minus1}(A)} p_Z(z)\ d\lambda^n(z)
%\notag \\
= \int_A p_f(x)\ d\lambda^n(x) \label{Background:CVF_classic}.
\end{align}

Since the equality (<ref>) is valid for all $\mathcal{A}_X$-measurable sets, the $\lambda^n$-unique probability density of $\mathbb{P}_f$ is given by

\begin{align}
p_f(x) = p_Z\left(f^{\minus1}(x)\right)\, \left|\det J_{f^{\minus1}}(x)\right|\quad (\text{a.e. } x\in X ).\label{Background:Density}
\end{align}

In order to unify the existing definitions of flows in the literature, we will speak of a flow or a proper flow when a density on $X$ of the form as in eq. (<ref>) is induced.

§ GENERALIZATION OF THE CVF

If it is argued that a flow indeed induces a probability density, often the CVF is merely named, and only in rare cases reference is made to sources like <cit.> or <cit.>. In most of the mathematical literature, it is proved in the following form

\begin{align}
\int_{f(Z)} \psi(x) dx = \int_Z \psi(f(z)) |\det J_f(z)| dz, \label{CVF_in_Literature}
\end{align}

where $f:U \to \R^n$ is injective, continuously differentiable with $U\subseteq \R^n$ open, $Z \subset U$ measurable, and $\psi$ Lebesgue integrable. Moreover, there are even broader formulations of this statement where $f$ is only differentiable or even only Lipschitz continuous everywhere and injective almost everywhere (cf. ). In <cit.> and <cit.>, generalizations are discussed where injectivity is not even required by considering the cardinality of the preimage set.

The identity (<ref>) is, however, rather analytically motivated for solving integrals and does not aim to provide a representation of the density induced by $f$. In the common case where $f$ satisfies Assumption <ref>, these two variants (eq. (<ref>) and (<ref>)) are valid, since both $f$ and $f^{\minus1}$ form diffeomorphisms. Nevertheless, when we consider generalizations, it is usually no longer clear how and whether they can also be applied for the purposes of normalizing flows. In the following sections, we derive a similar strong generalization of the CVF as in (<ref>), which is more suited to transforming probability densities, i.e., more suited for normalizing flows.

§.§ $\mathcal{L}$-Diffeomorphism

The basic idea is to require the conditions for a diffeomorphism only almost everywhere, since these sets do not affect integration. This idea leads to the following definition of a generalized transformation, providing a weaker set of conditions than in Assumption <ref>.
A mapping ${f:Z \to X}$ between two open sets $Z,X \subseteq \R^n$ is called Lebesgue-diffeomorphism ($\L$-diffeomorphism for short), if there are $\lambda^n$-null sets $N_Z, N_X$ with $N_Z$ closed such that the restriction $f:Z \setminus N_Z \to X \setminus N_X$ is bijective, continuously differentiable, and the set of critical points is a null set.

In the following, we list a few $\L$-diffeomorphism and their corresponding null sets:

* A diffeomorphism $f:Z \to X$ forms an $\L$-diffeomorphism where both $N_Z$ and $N_X$ are empty sets.

* The cubic function $f:\R \to \R$ with $f(x)=x^3$ is an $\L$-diffeomorphism with $N_Z=N_X=\emptyset$ and one critical point $\{0\}$, which is a null set.

* In particular, the plane polar coordinates transformation $f:\R_{+} \times [0,2\pi] \to \R^2$ given by $f(r,\phi) = (r \cos(\phi), r\sin(\phi))$ is an $\L$-diffeomorphism with

\begin{align*}
N_Z&=\{0\} \times [0,2\pi] \cup \R_{>0} \times \{0,2\pi\} \quad \text{and}\\
N_X&= \R_{+} \times \{0\}.
\end{align*}

In appendix Lemma <ref> we show that a$\L$-diffeomorphism is a measurable mapping with respect to the corresponding trace $\sigma$-algebras of the Lebesgue algebra, thus inducing a distribution on $X$ by the pushforward measure. Moreover, the name 'diffeomorphism' does justice to the Definition <ref>, which we state in the following lemma (see Appendix <ref> for the proof):

Let $f:Z \to X$ be an $\L$-diffeomorphism. Then there are $\lambda^n$-null sets $N_Z, N_X$ with $N_Z$ closed such that the restriction $f: Z \setminus N_Z \to X \setminus N_X$ is a diffeomorphism.

The reasoning behind the assumptions of an $\L$-diffeomorphism can be understood as follows: We can remove negligible sets from the domain and the target space such that the restriction is bijective and continuously differentiable. To show that the restricted inverse is continuously differentiable, we use the inverse function theorem (). Nevertheless, for the inverse function theorem to apply, all critical points and their image must still be removable, i.e., measure-theoretically negligible. That the set of critical points is a null set was assumed, and for the image, we use Sard's theorem (see Appendix <ref> for the detailed proof).

Furthermore, it is necessary to suppose that the set of critical points is a null set; because there are examples of continuously differentiable and bijective functions whose set of critical points does not have measure zero. These points cause the inverse function to be non-continuously differentiable on a set with a positive measure (see Appendix <ref> for an example).

§.§ CVF for $\L$-Diffeomorphism

The previously defined $\L$-diffeomorphisms form a reasonable generalization of diffeomorphisms and are comparable to those transformations discussed in the introduction of this section. Furthermore, the following theorem justifies the validity of CVF as well for $\L$-diffeomorphisms (see Appendix <ref> for the proof):
[CVF for $\L$-Diffeomorphism]
Let $f:Z \to X$ be an $\mathcal{L}$-diffeomorphism and $\mathbb{P}_Z$ a distribution on $Z$ with probability density $p_Z$ w.r.t. $\lambda^n$. Then the
CVF (see eq. (<ref>)) holds for $f$. In particular $f$ induces a distribution on $X$ with density $p_f$ given by

\begin{align}
p_f(x) = p_Z(f^{\minus1}(x))\, |\det J_{f^{\minus1}} (x)| &&(\text{a.e. } x \in X).
\end{align}

This theorem legitimizes the use of functions as proper flows that are not everywhere bijective, continuous, differentiable, or continuously differentiable. Even more, the inverse does not have to fulfill these properties everywhere either. In short, we can apply $\L$-diffeomorphisms as flows.

§.§ Invariance under Composition

The strength and tremendous upswing of normalizing flows mainly occurred because simple flows can be chained together. Thus, we can achieve the desired degree of complexity and expressiveness by increasing the number of simple flows. For this purpose, flows are often considered on the same base and target space. This crucial property is also retained for $\L$-diffeomorphisms (see Appendix <ref> for the proof):
Let $\Omega \subseteq \R^n$ be an open set and $f_1,f_2:\Omega \to \Omega$ $\mathcal{L}$-diffeomorphisms. Then the composition $f_2 \circ f_1$ is also an $\L$-diffeomorphism on $\Omega$.

From this lemma, it results inductively that the concatenation $f = f_K \circ \ldots \circ f_1$
of $K$ $\L$-diffeomorphisms $f_1,\ldots f_K$ on $\Omega$ forms an $\L$-diffeomorphism. Hence, common formulas from the normalizing flow literature, e.g., presented in <cit.> or <cit.>, also apply to $\L$-diffeomorphisms or can be extended to them. For example, the following holds for the Jacobian-determinant of $f$ with $x_i:= f_{i+1}^{\minus1} \circ \ldots \circ f_K^{\minus1}(x)$ and $x_K := x$

\begin{align*}
\det J_{f^{\minus1}}(x) = \prod_{i=1}^{K} \det J_{f_i^{\minus1}}(x_i) && (\text{a.e. } x \in \Omega).
\end{align*}

Despite the legitimate use of $\L$-diffeomorphisms mathematically, it is essential to note that these can lead to numerical instabilities. In some points, the Jacobian-determinant or the inverse function does not need to exist. Nevertheless, these values can be set meaningfully or ignored in some situations.

§ NON-SMOOTH ACTIVATIONS IN RESIDUAL FLOWS

In this section, we apply the previous results to residual mappings, which are perturbations of the identity of the form $f(x) = x + g(x)$. This justifies the use of non-smooth activations in planar, radial, and contractive residual flows.

§.§ Planar Flows

The term planar flow was first introduced by <cit.>, which refers to functions of the form

\begin{align*}
f_\text{P}(x) = x + u h(w^T x + b)
\end{align*}

with non-linearity $h$ and $w,u \in \R^n, b \in \R$. They describe a plane-wise expansion or contraction of all hyperplanes orthogonal to $w$. In order to admit also non-smooth activations, we generalize the results from <cit.> in the following theorem and obtain a sufficient criterion for the bijectivity of a planar flow $f_\text{P}$ (see Appendix <ref> for the proof).

Let $f_\text{P}$ be a planar flow with activation $h$. If the one-dimensional mapping $\psi:\R \to \R$ with

\begin{align*}
\psi(\lambda) = \lambda + w^T uh(\lambda)
\end{align*}

is bijective, then $f_\text{P}$ is also bijective.

However, bijectivity is not sufficient for the planar flow to induce a density of the desired form. But similar to Theorem <ref>, conditions can be imposed on a one-dimensional mapping such that $f_\text{P}$ describes an $\L$-diffeomorphism (see Appendix <ref> for the proof):

Let $f_\text{P}$ be a bijective planar flow. If there is a countable and closed set $N \subset \R$ such that the activation function $h$ is continuously differentiable on $\R \setminus N$ and

\begin{align*}
C= \left\{x \in \R \setminus N \mid 1 + w^Tu h'(x) = 0 \right\}
\end{align*}

is countable, then $f_\text{P}$ is an $\L$-diffeomorphism. In particular, $f_\text{P}$ is a proper flow.

Using Theorem <ref>, conditions for different activations such that the resulting planar flow describes a flow can be found via reducing it to a one-dimensional problem. The constraints for the most popular non-linearities are summarized in Table <ref>.

Conditions on the parameters of a planar flow $f_\text{P}$, such that it is a proper flow for the activation $h$ (see <ref> for the proofs).

Non-Linearity     Condition

ReLU         $w^Tu > -1$

ELU ($\alpha > 0$)      $w^Tu > \max(-1,-\tfrac{1}{\alpha})$

Tanh     $w^Tu \geq -1$

Softplus     $w^Tu > -1$

§.§ Radial Flows

Another intuitive way to perturb the identity in $\R^n$ is to expand or contract spherically around a centering point. This type of transformation was initially studied by <cit.> and subsequently by <cit.>. These transformations of the form

\begin{align*}
f_\text{R}(x) = x + \beta h\big(\|x-x_0\|_2\big) (x-x_0)
\end{align*}

are called radial flows, where $h:\R_+ \to \R_+$ is a localization function, $x_0 \in \R^n$ the center, and $\beta \in \R$. When $\beta$ is negative, a contraction occurs, and positive values lead to an expansion around the center $x_0$. The following theorem provides a sufficient criterion for the bijectivity of a radial flow if we consider non-smooth functions $h$ (see Appendix <ref> for the proof):

Let $f_\text{R}$ be a radial flow with localization $h$. If the one-dimensional mapping $\psi:\R_+ \to \R_+$ with

\begin{align*}
\psi(r) :=r + \beta h(r)r
\end{align*}

is bijective, then $f_\text{R}$ is also bijective.

Again, bijectivity is sufficient only for the existence of the inverse which means that $f_\text{R}$, in general, does neither describe an $\L$-diffeomorphism nor a flow. However, this property is ensured by the conditions of the following theorem (see Appendix <ref> for the proof):

Let $f_\text{R}$ be a bijective radial flow. If there is a countable, closed set $N \subset \R_{>0}$ such that the localization function $h$ is continuously differentiable on $\R_{>0} \setminus N$ and

\begin{align*}
C:= \left\{ r \in \R_{>0} \setminus N \left| 1+ \beta \left( h(r) + r h'(r)=0 \right) \right. \right\}
\end{align*}

is countable, then $f_\text{R}$ is an $\L$-diffeomorphism. In particular, $f_\text{R}$ is a proper flow.

§.§ Contractive Residual Flows

Contractive mappings provide a more general type of perturbations of the identity. A function $g$ is called contractive if there exists a constant $L<1$ such that for all $x,y \in \R^n$ it holds for any norm on the vector space $\R^n$

\begin{align*}
\|g(x) - g(y)\| \leq L\ \|x-y\|.
\end{align*}

In <cit.> and <cit.>, these kinds of residual flows are called (contractive) residual flows and are denoted henceforth by $f_\text{C}$. On the one hand, this strong condition on $g$ gives the bijectivity of $f_\text{C}$ by the Banach's fixed point theorem; on the other hand, it follows that $f_\text{C}$ has no critical points <cit.>. Thus, only a few assumptions are required to guarantee that a residual flow forms an $\L$-diffeomorphism.

Let $f_\text{C}$ be a residual flow with contractive perturbation $g$. If there is a closed $\lambda^n$-null set $N \subset \R^n$ such that $g$ is continuously differentiable on $\R^n \setminus N$, then $f_\text{C}$ is an $\L$-diffeomorphism. In particular, $f_\text{C}$ is a proper flow.

The proof of this statement follows directly from the property that Lipschitz continuous functions map $\lambda^n$-null sets to sets of $\lambda^n$-measure zero <cit.> and the fact that $\text{Lip}(f_\text{C}) = 1 + L$.

§ CONCLUSION

In this paper, we have shown that the conditions on a flow need not be strictly satisfied everywhere, and a certain degree of freedom is permitted on measure-theoretically negligible sets. With this, we have justified using $\L$-diffeomorphism instead of a normal diffeomorphism as flows. This gain significantly increases the possibilities in the design of normalizing flows and, in particular, allows the usage of non-smooth activations. Nevertheless, we have only justified their existence and usage mathematically, far from leading to a successful practical application. Thus, future work should investigate whether and in which situations non-smooth activations provide an actual gain. Moreover, we only applied these generalizations to simple residual flows. For this reason, it remains open to what extent other flows such as more general planar flows, like Sylvester flows <cit.>, or autoregressive flows <cit.> also benefit from this.

§ ACKNOWLEDGEMENTS

Niklas Koenen and Marvin N. Wright gratefully acknowledge the funding from the German Research Foundation (DFG) in the context of the Emmy Noether Grant 437611051.

§ GENERALIZATION OF THE CHANGE OF VARIABLES FORMULA

Let $f:Z \to X$ be an $\L$-diffeomorphism from the measurable space $(Z,\mathcal{A}_Z)$ into the target space $(X, \mathcal{A}_X)$. Then $f$ is an $\mathcal{A}_Z$-$\mathcal{A}_X$-measurable mapping, i.e., a random variable.

According to the definition of measurable functions, we have to show that every preimage of a measurable set contained in $\mathcal{A}_X$ is also an element of $\mathcal{A}_Z$. Let $B \in \mathcal{A}_X$. Since $f$ is an $\L$-diffeomorphism, there exist $\lambda^n$-null sets $N_Z, N_X$ such that the restriction $f:Z \setminus N_Z \to X \setminus N_X$ is bijective and continuously differentiable. We can decompose the preimage of $B$ as

\begin{align*}
f^{\minus1}(B) = \left(N_Z \cap f^{\minus1}(B)\right) \cup \left(Z \setminus N_Z \cap f^{\minus1}(B)\right).
\end{align*}

Since $N_Z$ is a null set, the subset $N_Z \cap f^{\minus1}(B)$ also has a measure of zero. We can represent the other set of the equation above by the bijectivity of $f$ on $Z \setminus N_Z$ as

\begin{align}\label{Appendix:Lemma_1_hard_set}
Z \setminus N_Z \cap f^{\minus1}(B) = f^{\minus1}(X \setminus N_X \cap B).
\end{align}

Thus it results from the continuity of the restricted mapping, which is in particular $\L(Z \setminus N_Z)$-$\L(X\setminus N_X)$-measurable, that the set (<ref>) is contained in the trace $\sigma$-algebra ${\L(Z \setminus N_Z)}$, hence an element of $\mathcal{A}_Z$. In the end, we can represent $f^{\minus1}(B)$  as a union of a null set and a $\mathcal{A}_Z$-measurable set, which shows the claim.

From the Definition <ref> of an $\L$-diffeomorphism, we obtain the existence of null sets $N_Z, N_X$ with $N_Z$ closed such that the restricted function $f:Z \setminus N_Z \to X \setminus N_X$ is bijective, continuously differentiable, and the set of critical points $C$ has measure zero. At first, we show that $C$ is a closed set. Consider the mapping $g:Z \setminus N_Z \to \R$ with $g(x) = \det(J_f(x))$. Because of the continuity of the determinant and the continuous differentiability of the restricted $f$, the mapping $g$ is also continuous. Therefore the preimage under $g$ of closed sets is also closed according to the topological definition of continuity, i.e.,

\begin{align*}
g^{\minus1}(\{0\}) = \{x \in Z \setminus N_Z \mid \det J_f(x)=0\} = C
\end{align*}

is a closed set; thus, $Z \setminus (N_Z \cup C)$ is open in $\R^n$. It follows by the inverse function theorem and its consequences <cit.> that the restriction

\begin{align*}
f:Z \setminus (N_Z \cup C) \to X \setminus (N_X \cup f(C))
\end{align*}

forms a diffeomorphism. By assumption, $N_C$ and $C$ are $\lambda^n$-null sets, and so $N_Z \cup C$. From Sard's theorem <cit.>, it follows that $f(C)$, hence especially $N_X \cup f(C)$, is a set of measure zero.

We obtain from Lemma <ref> the existence of two null sets $N_Z, N_X$ so that the restriction $f:Z \setminus N_Z \to X \setminus N_X$ of an $\L$-diffeomorphism  forms a diffeomorphism. Let $A \in \mathcal{A}_X$ be a measurable set in the target space $X$. Because of the additivity of $\mathbb{P}_Z$ and the definition of the pushforwad measure $\mathbb{P}_f$, it follows

\begin{align}
\mathbb{P}_f(B):= \mathbb{P}_Z\left(f^{\minus1}(B)\right)
&= \mathbb{P}_Z\left((f^{\minus1}(B)\setminus N_Z)\ \dot\cup\ (f^{\minus1}(B)\cap N_Z) \right)\notag \\
&= \mathbb{P}_Z\left(f^{\minus1}(B)\setminus N_Z\right) + \mathbb{P}_Z\left(f^{\minus1}(B)\cap N_Z \right). \label{Appendix:Theorem_CVF_sets}
\end{align}

Since $\mathbb{P}_Z$ is absolutely continuous w.r.t. the Lebesgue measure and $f^{\minus1}(B) \cap N_Z$ is a subset of a $\lambda^n$-null set, the right summand of term (<ref>) vanishes. For the other term, we note the following set equality:

\begin{align}
f^{\minus1}(B) \setminus N_Z = f^{\minus1} (B \setminus N_X).
\end{align}

Because the restriction $f:Z\setminus N_Z \to X \setminus N_X$ is a diffeomorphism and $B\setminus N_X$ is an element of the trace $\sigma$-algebra $\L(X \setminus N_X)$, the equation (<ref>) can be extended by the CVF from eq. (<ref>) as

\begin{align}
\mathbb{P}_f(B) &= \int_{f^{\minus1}(B \setminus N_X)} p_z(z)\ d\lambda^n(z)
= \int_{B \setminus N_X} p_z\left(f^{\minus1}(x)\right)\ \left|\det J_{f^{\minus1}}(x)\right|\ d\lambda^n(x). \label{Appendix:Theorem_CVF_integral}
\end{align}

We define a function $p_f: X \to \R_+$ with

\begin{align*}
p_f(x) = \left\{
\begin{array}{cl}
p_Z\left(f^{\minus1}(x)\right)\ \left|\det J_{f^{\minus1}}(x)\right| & x \in X \setminus N_X \\
\vphantom{\Big|}0 & x \in N_X
\end{array} \right..
\end{align*}

Since $B \cap N_X \subseteq N_X$ is a $\lambda^n$-null set, we can additionally integrate in (<ref>) over this set and obtain the expression of the distribution induced by $f$

\begin{align}
\mathbb{P}_f(B)= \int_{B \setminus N_X} p_f(x)\ d\lambda^n(x)  + \int_{B \cap N_X} p_f(x)\ d\lambda^n(x)= \int_{B} p_f(x)\ d\lambda^n(x),
\end{align}

hence $p_f$ forms a density of distribution $\mathbb{P}_f$. Finally, by the Radon-Nikodym theorem <cit.>, the density function is $\lambda^n$-unique, which shows the claim.

Since $f_1$ and $f_2$ are $\L$-diffeomorphisms, by Lemma <ref> there are $\lambda^n$-null sets $N_Z^1, N_Z^2, N_X^1, N_X^2 \subset \Omega$ with $N_Z^1$ and $N_Z^2$ closed such that both

\begin{align}
f_1 : \Omega \setminus N_Z^1 \to \Omega \setminus N_X^1\quad \quad \text{ and } \quad \quad
f_2 : \Omega \setminus N_Z^2 \to \Omega \setminus N_X^2 \label{Appendix:Lemma_Invariance_diffeos}
\end{align}

are diffeomorphisms. We remove from the domain of $f_1$ all elements mapping to the set $N_Z^2$ since $f_2$ is not a diffeomorphism on it. Hence we define

\begin{align}
N_Z := N_Z^1 \cup f_1^{\minus1} \left(N_Z^2 \cap \Omega \setminus N_X^1\right) \quad \quad \text{ and } \quad \quad
N_X := N_X^2 \cup f_2\left(N_X^1 \cap \Omega \setminus N_Z^2\right). \label{Appendix:Lemma_Invariance_1}
\end{align}

It is relatively easy to see that both $f_1(\Omega \setminus N_Z) = \Omega \setminus (N_X^1 \cup N_Z^2)$ and $f_2(f_1(\Omega \setminus N_Z)) = \Omega \setminus N_X$ are valid. Moreover, from the continuity of $f_1$ and since $N_Z^2 \cap \Omega \setminus N_X^1$ is closed in the subspace topology on $\Omega \setminus N_X^1$, it follows that $f^{\minus1}(N_Z^2 \cap \Omega \setminus N_X^1)$ is closed in $\Omega \setminus N_Z^1$. Hence there exists a closed set $A \subset \Omega$ such that (see Sec. 6 in  for more information)

\begin{align*}
f^{\minus1}\left(N_Z^2 \cap \Omega \setminus N_X^1\right)= A \cap \Omega \setminus N_Z^1.
\end{align*}

This results with eq. (<ref>) in

\begin{align*}
N_Z = N_Z^1 \cup \left(A \cap \Omega \setminus N_Z^1\right) = N_Z^1 \cup A
\end{align*}

which is a closed set in $\Omega$. In summary, the mapping $f_2 \circ f_1 : \Omega \setminus N_Z \to \Omega \setminus N_X$ is well-defined, continuously differentiable, bijective, and has no critical points, since we have merely further restricted the diffeomorphisms from (<ref>). In addition, $N_Z$ is a closed set in $\Omega$. To complete the proof, we still need to show that $N_Z$ and $N_X$ have a $\lambda^n$-measure of zero. However, this follows directly from the fact that diffeomorphisms map null sets to null sets.

The following example illustrates why we assume that the set of critical points is a null set. Indeed, one can construct bijective and continuously differentiable functions whose inverse is not differentiable on a set with positive measure.

As an example, consider the following recursively defined set on the interval $[0,1]$ called the Smith-Volterra-Cantor set $C_\text{SV}$. The definition and properties of this set are only briefly sketched here. For more information and detailed proofs, the reader is referred to <cit.> or <cit.>. Another example can also be found in <cit.>.

The recursive definition of $C_\text{SV}$ on $[0,1]$ starts by removing the open middle quarter from the interval. Then we create subsequent sets $S_n$ by removing an open interval of length $\tfrac{1}{4^n}$ from the center of each interval in $S_{n\minus1}$, i.e.,

$S_1=$        $\Big[0, \frac{3}{8} \Big]$        $\cup$        $\Big[\frac{5}{8}, 1 \Big]$

$S_2=$    ${\small \Big[0, \frac{5}{32} \Big]}$     $\cup$     ${\small \Big[\frac{7}{32}, \frac{3}{8} \Big]}$     $\cup$     ${\small \Big[\frac{5}{8}, \frac{25}{32} \Big]}$     $\cup$     ${\small \Big[\frac{27}{32}, 1 \Big]}$

$S_3=$    $ {\tiny \Big[0, \tfrac{1}{16} \Big] \cup \Big[\tfrac{3}{32}, \tfrac{5}{32} \Big]}$     $\cup$     ${ \tiny \Big[\tfrac{7}{32}, \frac{9}{32} \Big] \cup \Big[\tfrac{5}{16}, \frac{3}{8} \Big]}$    $\cup$
$ {\tiny \Big[\tfrac{5}{8}, \tfrac{11}{16} \Big] \cup \Big[\tfrac{23}{32}, \tfrac{25}{32} \Big]}$    $\cup$
${\tiny \Big[\tfrac{27}{32}, \tfrac{29}{32} \Big] \cup \Big[\tfrac{15}{16}, 1 \Big]}$

[decoration=Cantor set,line width=3mm]
(0,0) – (1,0);

(0,-0.5) – (3/8, -0.5);
(5/8,-0.5) – (1, -0.5);

(0,-1) – (5/32, -1);
(7/32,-1) – (3/8, -1);
(5/8,-1) – (25/32, -1);
(27/32,-1) – (1, -1);

(0,-1.5) – (2/32, -1.5);
(3/32,-1.5) – (5/32, -1.5);
(7/32,-1.5) – (9/32, -1.5);
(10/32,-1.5) – (12/32, -1.5);
(20/32,-1.5) – (22/32, -1.5);
(23/32,-1.5) – (25/32, -1.5);
(27/32,-1.5) – (29/32, -1.5);
(30/32,-1.5) – (1, -1.5);

(0,-2) – (2/64 - 0.25^4, -2);
(2/64 + 0.25^4,-2) – (2/32, -2);
(3/32,-2) – (4/32 - 0.25^4, -2);
(4/32 + 0.25^4,-2) – (5/32, -2);
(7/32,-2) – (8/32 - 0.25^4, -2);
(8/32 + 0.25^4,-2) – (9/32, -2);
(10/32,-2) – (11/32 - 0.25^4, -2);
(11/32 + 0.25^4,-2) – (12/32, -2);
(20/32,-2) – (21/32 - 0.25^4, -2);
(21/32 + 0.25^4,-2) – (22/32, -2);
(23/32,-2) – (24/32 - 0.25^4, -2);
(24/32 + 0.25^4,-2) – (25/32, -2);
(27/32,-2) – (28/32 - 0.25^4, -2);
(28/32 + 0.25^4,-2) – (29/32, -2);
(30/32,-2) – (31/32 - 0.25^4, -2);
(31/32 + 0.25^4,-2) – (1, -2);

figureVisualization of $S_n$.

Now the Smith-Volterra-Cantor set is defined as the following intersection of all $S_n$

\begin{align*}
C_\text{SV} := \bigcap_{n=1}^\infty S_n,
\end{align*}

which is closed because of the countable intersection of closed sets. In addition, in each recursion step in each of the $2^{n\minus1}$ intervals, the middle pieces with the Lebesgue measure of $\tfrac{1}{4^n}$ are removed. Hence it holds

\begin{align*}
\lambda^1(C_\text{SV}) = 1 - \sum_{n=1}^\infty \frac{2^{n\minus1}}{4^n} = 1 - \frac{1}{4} \sum_{n=0}^\infty \left(\frac{1}{2} \right)^n =1 - \frac{1}{2} = \frac{1}{2}
\end{align*}

due to the geometric series. In order to construct a counterexample, we consider the function $d:[0,1] \to \R$ with $d(x)= \inf_{c \in C_\text{SV}} |x-c|$. This function is obviously continuous and $d(x) = 0$ holds if and only if $x \in C_\text{SV}$, otherwise it only takes positive values. Because of the continuity, this function is integrable and it follows from the fundamental theorem of calculus that

\begin{align*}
f:(0,1) &\to f\big((0,1)\big) \\
f(x) &= \int_0^x d(z)\ dz,
\end{align*}

is continuously differentiable. Moreover, $f$ is injective resulting from the construction of $C_\text{SV}$ and the fact that $d$ is not constant zero on any interval with positive measure. But the set of critical points of $f$ is the Smith-Volterra-Cantor set, which is not a $\lambda$-null set.

§ NON-SMOOTH ACTIVATIONS IN RESIDUAL FLOWS

§.§ Planar Flows

Theorem <ref> presented here corresponds to a generalization of the proof given by <cit.> since not only the smooth activation function $h(x)=\tanh(x)$ is considered. However, the argument in the proof is very similar.

Let $y \in \R^n$ be arbitrary. Then we have to show that the following equation has a unique solution

\begin{align}
f_\text{P}(x) = x + u h\left(w^T x + b\right) = y. \label{Appendix:Planar_Flows_original_eq}
\end{align}

If $w=0$ holds, then a unique solution is given by $x = y - uh(b)$, wherefrom bijectivity results. For this reason, let $w \neq 0$ in the following, thus $w$ spans a one-dimensional linear subspace $W$ in $\R^n$. Consequently, each element $x \in \R^n$ has a unique orthogonal decomposition $x = x_\parallel + x_\perp$ with $x_\parallel \in W$ and $x_\perp \in W^\perp := \{x \in \R^n \mid w^Tx = 0\}$. Due to the orthogonality of $x_\perp$ and $w$, the following solution of the orthogonal component depending on the parallel one can be inferred from eq. (<ref>)

\begin{align}
x_\perp = y - x_\parallel - uh\left(w^Tx_\parallel + b\right).
\end{align}

Since $W$ is a one-dimensional linear subspace, there is a unique $\lambda \in \R$ with $x_\parallel = w \tfrac{\lambda}{w^Tw}$. By using this representation, the original equation multiplied by $w^T$ from the left yields

\begin{align}
w^T y &= w^T x_\perp + w^Tw \frac{\lambda}{w^Tw} + w^T u h(\lambda + b) = \lambda + w^T u h(\lambda + b) \label{Appendix:Planar_Flows_eq_1}\\
&= \psi(\lambda + b) - b, \notag
\end{align}

where the last equation in (<ref>) follows again by $w^Tx_\perp = 0$. The assumed bijectivity of $\psi$ leads to the unique existence of a $\lambda_y \in \R$, which solves the equation above. In addition, this implies the existences of $x_\parallel$ and $x_\perp$, thus

\begin{align*}
x = x_\parallel + x_\perp = y - uh\left(\lambda_y + b\right)
\end{align*}

is the unique solution of the initial equation (<ref>) and consequently the bijectivity of $f_\text{P}$ follows.

In case $w=0$, $f_\text{P}$ represents an affine linear mapping describing a diffeomorphism, hence an $\L$-diffeomorphism. For this reason, let $w \neq 0$. Consequently, $w$ spans a one-dimensional linear subspace $W$, from which follows a unique orthogonal decomposition of the vector space $\R^n = W + W^\perp$ as a direct sum. This leads to a characterization of the vector space as a disjoint union of hyperplanes, i.e.,

\begin{align*}
\R^n = \dot{\bigcup_{\lambda \in  \R}} H_\lambda \quad\quad \text{with} \quad\quad H_\lambda := \left\{ \left. w \frac{\lambda}{w^Tw} + x_\perp \right| x_\perp \in W^\perp \right\}.
\end{align*}

Consider the function $\tau : \R^n \to \R$ with $\tau(x) = w^Tx + b$ mapping every element of a hyperplane $H_\lambda$ to the same value $\lambda + b$. Since the activation function $h$ is not continuously differentiable on the countable set $N$, we remove all hyperplanes mapping to $N$ under $\tau$. So we define

\begin{align*}
H_N := \bigcup_{n \in N} H_ {n - b} \quad\quad \text{such that} \quad\quad \tau(H_N) = N.
\end{align*}

Because $\tau$ is continuous and $N$ is closed in $\R$, it follows that $H_N = \tau^{\minus1}(N)$ is closed in $\R^n$. Moreover, $H_N$ as a countable union of hyperplanes is also a null set w.r.t. the Lebesgue measure $\lambda^n$ due to the subadditivity. Accordingly, the restriction $f_\text{P}:\R^n \setminus H_N \to \R^n \setminus f_\text{P}(H_N)$ is bijective and continuously differentiable. Furthermore, $H_N$ is a closed $\lambda^n$-null set. We note that the image of a hyperplane $H_\lambda$ under $f_\text{P}$ is also a hyperplane; hence $f_\text{P}(H_N)$ forms a countable union of hyperplanes which is a null set. Finally, it is left to show that the set of critical points of the restricted planar flow is also a null set. By the matrix determinant lemma we get for the set of critical points the set equality

\begin{align*}
C_{f_\text{P}}:=\left\{x \in \R^n \setminus H_N \left| \det J_{f_\text{P}}(x) = 0 \right.\right\}
= \left\{x \in \R^n \setminus H_N \left| 1 + w^Tu h'(\tau(x)) = 0 \right.\right\}.
\end{align*}

This gives $\tau(C_{f_\text{P}})=C$, which we have assumed to be a countable set. Hence

\begin{align*}
C_{f_\text{P}} = \tau^{\minus1}(C)= \bigcup_{c \in C} \tau^{\minus1}(\{c\}) = \bigcup_{c \in C} H_{c - b},
\end{align*}

which is as a countable union of hyperplanes a $\lambda^n$-null set.

§.§.§ Examples for some Activations

In the following, we infer conditions for particular choices of activation functions for a planar flow $f_\text{P}$. A visualization of the crucial function for each proof can be found in Figure <ref>, where the conditions on the flow are fulfilled, just (not) fulfilled and not fulfilled anymore.

If $w^Tu>-1$ holds, then $f_\text{P}$ with activation $\text{ReLU}$ is an $\L$-diffeomorphism.

Consider the mapping $\psi:\R \to \R$ given by

\begin{align*}
\psi(\lambda) = \lambda + w^Tu\, \text{ReLU}(\lambda + b) = \left\{\begin{array}{ll}
\lambda &, \lambda < -b  \\
\lambda + w^Tu (\lambda + b)     &, \lambda \geq -b
\end{array} \right..
\end{align*}

Because of the inequality $w^Tu > -1$, this mapping is strictly monotonically increasing and thus obviously bijective. Consequently, the bijectivity of the planar flow $f_\text{P}$ follows from Theorem <ref>. Furthermore, the $\text{ReLU}$ activation is continuously differentiable on $\R \setminus \{0\}$ and it holds for all $\lambda \in \R \setminus \{0\}$

\begin{align*}
1 + w^Tu\, \text{ReLU}'(\lambda)= \left\{ \begin{array}{ll}
1 &, \lambda < 0  \\
1 + w^Tu &, \lambda > 0
\end{array} \right. \neq 0.
\end{align*}

Therefore, the planar flow $f_\text{P}$ with activation $\text{ReLU}$ has no critical points; thus, the claim follows from Theorem <ref>.

If $w^Tu\geq-1$ holds, then $f_\text{P}$ with activation $\tanh$ is an $\L$-diffeomorphism.

Consider the function $\psi: \R \to \R$ given by $\psi(\lambda) =  \lambda + w^Tu \tanh(\lambda)$ with derivative $\psi'(\lambda)= 1 + \frac{w^Tu}{\cosh(\lambda)^2}$. Since the hyperbolic cosine is $1$ only at $0$ and otherwise always greater than $1$, we get with the assumed inequality $\psi'(\lambda) \geq 0$ and equality only if $\lambda=0$ and $w^Tu = -1$. Therefore, the function $\psi$ is strictly monotonically increasing, hence injective. Moreover, the surjectivity follows from the boundedness of the hyperbolic tangent. Consequently, the bijectivity of the planar flow $f_\text{P}$ follows from Theorem <ref>. The activation function is continuously differentiable, and as seen earlier, the equation $1+w^Tu \tanh(\lambda)=0$ is only satisfied if $w^Tu=-1$ and $\lambda = 0$. In any case, the set of critical points is a countable set, so the claim follows from Theorem <ref>.

If $w^Tu> \max(-1, -\frac{1}{\alpha})$ holds, then $f_\text{P}$ with activation $\text{ELU}$ $(\alpha > 0)$ is an $\L$-diffeomorphism.

Consider the function $\psi:\R \to \R$ with

\begin{align*}
\psi(\lambda) = \lambda + w^Tu\, \text{ELU}(\lambda +b) = \left\{ \begin{array}{ll}
\lambda + w^Tu \alpha \left( e^{\lambda+b} -1 \right)     &, \lambda \leq -b  \\
\lambda + w^Tu(\lambda + b)     &, \lambda > -b
\end{array} \right..
\end{align*}

This function is continuously differentiable on $\R \setminus \{b\}$ with derivative given by

\begin{align*}
\psi'(\lambda) = \left\{ \begin{array}{ll}
1 + w^Tu \alpha e^{\lambda+b}     &, \lambda < -b  \\
1 + w^Tu     &, \lambda > -b
\end{array} \right..
\end{align*}

By the assumed inequality and $e^{\lambda + b} \in (0,1)$ for $\lambda < -b$, we obtain for $\lambda < -b$

\begin{align}
\psi'(\lambda)= 1+ w^Tu \alpha e^{\lambda + b} > 1+ \max\left(-1,- \frac{1}{\alpha} \right)\alpha e^{\lambda + b} > 1 + \max(-\alpha, -1) \geq 0. \label{PlanarFlow:ELU_1}
\end{align}

In addition, we get the positivity for the other case $\lambda > -b$

\begin{align}
\psi'(\lambda) = 1 + w^Tu > 1+ \max\left( -1, -\frac{1}{\alpha} \right) \geq 0. \label{PlanarFlow:ELU_2}
\end{align}

Because of the limit $\lim_{\lambda \searrow -b} \psi(\lambda) = \psi(-b)= -b$, the function $\psi$ is strictly monotonically increasing, thus injrective. Furthermore, surjectivity results from the mean value theorem; hence $f_\text{P}$ is bijective by Theorem <ref>. Additionally, the activation function $\text{ELU}$ is continuously differentiable on $\R \setminus \{b\}$, and similar to (<ref>) and (<ref>) the inequality $\psi'(\lambda) = 1 + w^Tu \text{ELU}(\lambda)> 0$ follows for all $\lambda \neq 0$. Therefore we can apply Theorem <ref>, which shows the claim.

axis x line=center,
axis y line=center,
xlabel style=below,
ylabel style=left,
legend pos= north west,
legend cell align = left,
legend style=nodes=scale = 0.75, line width = 1]
[mark=none, color=green!50!black, line width=1 ,domain=-3.5:0] x;

[mark=none, color=blue!50!black, line width=1 ,domain=-3.5:0] x;

[mark=none, color=magenta!50!black, line width=1 ,domain=-3.5:0] x;
$w^Tu = 1$, $w^Tu = -1$, $w^Tu = -2$;

[mark=none,color=green!50!black, line width=1, domain=0:3.5, smooth] 2*x;
[mark=none,color=blue!50!black, line width=1, domain=0:3.5, smooth] 0;
[mark=none,color=magenta!50!black, line width=1, domain=0:3.5, smooth] -x;

$\text{ReLU}$ activation

samples = 200,
axis x line=center,
axis y line=center,
xlabel style=below,
ylabel style=left,
legend pos= north west,
legend cell align = left,
legend style=nodes=scale = 0.75, line width = 1]
[mark=none, color=green!50!black, line width=1 ,domain=-3.5:3.5] x + tanh(x);

[mark=none, color=blue!50!black, line width=1 ,domain=-3.5:3.5] x - tanh(x);

[mark=none, color=magenta!50!black, line width=1 ,domain=-3.5:3.5] x - 2* tanh(x);
$w^Tu = 1$, $w^Tu = -1$, $w^Tu = -2$;

$\tanh$ activation

axis x line=center,
axis y line=center,
xlabel style=below,
ylabel style=left,
legend pos= north west,
legend cell align = left,
legend style=nodes=scale = 0.75, line width = 1]
[mark=none, color=green!50!black, line width=1 ,domain=-3.5:0] x + 2*(exp(x) - 1);

[mark=none, color=blue!50!black, line width=1 ,domain=-3.5:0] x - 0.5* 2*(exp(x) - 1);

[mark=none, color=magenta!50!black, line width=1 ,domain=-3.5:0] x - 0.9* 2* (exp(x) - 1);
$w^Tu = 1$, $w^Tu = -0.5$, $w^Tu = -0.9$;

[mark=none,color=green!50!black, line width=1, domain=0:3.5, smooth] x + x;
[mark=none,color=blue!50!black, line width=1, domain=0:3.5, smooth] x - 0.5*x;
[mark=none,color=magenta!50!black, line width=1, domain=0:3.5, smooth] x - 0.9*x;

$\text{ELU}$ activation with $\alpha = 2$

samples = 200,
axis x line=center,
axis y line=center,
xlabel style=below,
ylabel style=left,
legend pos= north west,
legend cell align = left,
legend style=nodes=scale = 0.75, line width = 1]
[mark=none, color=green!50!black, line width=1 ,domain=-3.5:3.5] x + 0.5*ln(1 + exp(x));

[mark=none, color=blue!50!black, line width=1 ,domain=-3.5:3.5] x - ln(1+exp(x));

[mark=none, color=magenta!50!black, line width=1 ,domain=-3.5:3.5] x - 2* ln(1 + exp(x));
$w^Tu = 0.5$, $w^Tu = -1$, $w^Tu = -2$;

$\text{Softplus}$ activation

Visualization of the function $\psi$ considered in the proofs of Section <ref>. For each of the activations $\text{ReLU}, \tanh, \text{ELU}(\alpha = 2)$ and $\text{Softplus}$, we plotted the function $\psi$ with $b=0$ for three different choices of $w^Tu$; the conditions for an $\L$-diffeomorphism are satisfied (green), just (not) satisfied (blue), and not satisfied (magenta).

If $w^Tu> -1$ holds, then $f_\text{P}$ with activation $\text{Softplus}$ is an $\L$-diffeomorphism.

For the function $\psi:\R \to \R$ given by $\psi(\lambda) =  \lambda + w^Tu \text{Softplus}(\lambda + b)$ with $\text{Softplus}(\lambda)= \log\left(1+e^{\lambda + b} \right)$ the derivative is

\begin{align*}
\psi'(\lambda) = 1 + w^Tu \frac{e^{\lambda + b}}{1 + e^{\lambda + b}} = 1 + w^Tu \frac{1}{1 + e^{\minus(\lambda + b)}}.
\end{align*}

Since the range of the factor $\frac{1}{1 + e^{\minus(\lambda + b)}}$ is in the interval $(0,1)$, the assumed inequality results for all $\lambda \in \R$ in

\begin{align}
\psi'(\lambda) > 1- \frac{1}{1 + e^{\minus(\lambda + b)}} > 0.
\end{align}

Consequently, the function $\psi$ is strictly monotonically increasing, hence injective. By continuity, the limit $\lim_{\lambda \to - \infty} = -\infty$ results. Without loss of generality, we consider $\lambda > -b$ in the following. In this case, the inequality $\log(1 + e^{\lambda +b}) \leq \log(2 e^{\lambda +b})$ holds, and we obtain for $w^Tu < 0$

\begin{align*}
\lim_{\lambda \to \infty} \psi(\lambda) \geq \lim_{\lambda \to \infty} \lambda + w^Tu \log\left(2e^{x+b}\right) = \lim_{\lambda \to \infty} \lambda + w^Tu \log(2) \left(\lambda + b\right) = \infty.
\end{align*}

For the other case $w^Tu \geq 0$, results $\lim_{\lambda \to \infty} \psi(\lambda) \geq \lim_{\lambda \to \infty} \lambda = \infty$. Thus, the surjectivity follows from the mean value theorem. Therefore, we get the bijectivity from Theorem <ref>. In addition, the activation function $\text{Softplus}$ is everywhere continuously differentiable and has not critical points, so the claim follows from Theorem <ref>.

§.§ Radial Flows

Theorem <ref> presented here corresponds to a generalization of the proof given by <cit.>, since not only the localization function $h(r)=\tfrac{1}{\alpha + r}$ is considered. However, the arguments in the proof are very similar.

Let $y \in \R^n$ arbitrary. We have to show for bijectivity that the following equation has a unique solution

\begin{align}
f_\text{R}(x) = x + \beta h(r) (x-x_0)=y \quad\quad \text{where} \quad\quad r:= \|x-x_0\|_2. \label{Appendix:Radial_Bijectivity_1}
\end{align}

If $y=x_0$, we then obtain the equality

\begin{align*}
0 = \left\|y-x_0\right\|_2 = \left\|x-x_0 + \beta h(r) (x-x_0) \right\|_2 = r + \beta h(r)r = \psi(r)
\end{align*}

after rearranging equation (<ref>) and taking the Euclidean norm. Since $\psi(0) = 0$ and $\psi$ was assumed to be bijective, it follows that $\|x-x_0\|=0$ which is equivalent to $x=x_0$ because of the definiteness of the Euclidean distance. Hereafter let $y \in \R^n \setminus \{x_0\}$, i.e., $r=\|x-x_0\|_2> 0$. In this case, each element $x - x_0$ with $x \in \R^n \setminus \{x_0\}$ can be expressed unambiguously as the product of its projection on the unit sphere and its Euclidean distance. Thus, for $x \in \R^n \setminus \{x_0\}$ there exists a unique normalized direction vector $\hat{x} \in S_1(x_0) := \{x \in \R^n \mid \|x-x_0\|_2 = 1\}$ such that $ x= x_0 + r\hat{x}$ where $r:= \|x-x_0\|_2$. Substituting this representation into equation (<ref>), one obtains after conversion

\begin{align}
y-x_0 = \hat{x} \left( r + \beta h(r) r \right) = \hat{x} \psi(r), \label{Appendix:Radial_Bijectivity_2}
\end{align}

and by taking of the norm of this

\begin{align*}
0 < \|y-x_0\|_2 = \|\hat{x} \psi(r)\|_2 = \psi(r).
\end{align*}

Because of the bijectivity of $\psi$, the unique existence of a radius $r_y > 0$ around the centering point results. Due to the fact that $\psi(r)>0$ remains valid for $r>0$, equation (<ref>) gives the following unambiguous expression of the direction vector $\hat{x}$

\begin{align*}
\hat{x}_y := \frac{y-x_0}{\psi(r_y)} = \frac{y-x_0}{r_y + \beta h(r_y) r_y}.
\end{align*}

Hence $x = x_0 + r_y \hat{x}_y$ is the unique solution of the original equation (<ref>), concluding finally that $f_\text{R}$ is bijective.

According to the requirement, the localization function is not continuously differentiable for all radii, so spheres with such distances around the centering point must be removed from the domain of $f_\text{R}$. Furthermore, the point $x_0$ corresponding to a radius of $0$ must be eliminated in order to restrict the localization function to an open set, thus allowing us to verily speak of differentiability. For this purpose we define

\begin{align*}
S:= \bigcup_{r \in N} S_r(x_0) \cup \{x_0\} \quad\quad \text{where}\quad\quad S_r(x_0):=\left\{ x \in \R^n \mid \|x-x_0\|_2=r  \right\}.
\end{align*}

Since the shifted Euclidean norm $\tau(x) := \|x-x_0\|_2 $ is continuous, it follows that the preimage of the closed set $N \cup \{0\}$ is also closed. Therefore, the set of eliminated points

\begin{align*}
S= \bigcup_{r \in N} \tau^{\minus1}\left(\{r\}\right) \cup \tau^{\minus1}(\{0\})= \tau^{\minus1}\left(N \cup \{0\}\right)
\end{align*}

is closed; moreover, it is a set of measure zero because a countable union of spheres and points is a $\lambda^n$-null set. In summary, the restriction $f_\text{R}: \R^n \setminus S \to \R^n \setminus f_\text{R}(S)$ describes a bijective and continuously differentiable mapping. Besides, the following is valid for every $x\in S_r(x_0)$

\begin{align}
\|f_\text{R}(x) - x_0 \|_2 = \left\| x - x_0 + \beta h(r) (x-x_0) \right\|_2 = r+ \beta h(r)r. \label{Appendix:Radial_Main_1}
\end{align}

This indicates that $f_\text{R}$ maps spheres with radius $r$ to spheres with radius $r+ \beta h(r)r$ around the centering point $x_0$; thus, $f_\text{R}(S)$ is also a countable union of null sets, therefore, itself a null set. Finally, it remains to show that the set of critical points of this restriction is also a Lebesgue null set. From the matrix determinant lemma and the higher dimensional differentiation rules, the Jacobian-determinant of $f_\text{R}$ at position $x\in \R^n \setminus S$ is given by

\begin{align}
\det\left(J_{f_\text{R}}(x)\right) = \left(1+\beta h(r)\right)^{n-1} \left(1+\beta h(r) + \beta h'(r) r\right). \label{Appendix:Radial_Main_det}
\end{align}

Since $f_\text{R}(x_0)=x_0$ holds and $f_\text{R}$ is bijective, both $\|f_\text{R}(x) - x_0 \|_2 > 0$ and $r = \|x - x_0\|_2>0$ hold for all $x \in \R^n \setminus S$. Thus equation (<ref>) yields

\begin{align*}
\frac{\|f_\text{R}(x) - x_0 \|_2}{r} = 1 + \beta h(r) > 0.
\end{align*}

Knowing that this term does not vanish for any $x\in \R^n \setminus S$, the set of critical points can be represented using equation (<ref>) like

\begin{align*}
C_{f_\text{R}} := \left\{x \in \R^n \setminus S \left| \det J_{f_\text{R}}(x) = 0 \right.\right\} = \left\{x \in \R^n \setminus S \left| 1 + \beta h(\tau(x))+ \beta h'(\tau(x)) = 0 \right.\right\}.
\end{align*}

This gives $\tau(C_{f_\text{R}}) = C$, which we have assumed to be a countable set. Hence

\begin{align*}
C_{f_\text{R}} = \tau^{\minus1}(C) = \bigcup_{c \in C} \tau^{\minus1}(\{c\}) = \bigcup_{c \in C} S_c(x_0),
\end{align*}

which is as a countable union of spheres a $\lambda^n$-null set.


# Deep Sequential Models for Suicidal Ideation from Multiple Source Data
††thanks:  Manuscript received December 17, 2018. This work has been funded in
part by the Spanish MINECO under Grants TEC2015-69868-C2-1-R,
TEC2016-78434-C3-3-R and TEC2017‐92552‐EXP, in part by Spanish MICINN under
Grant RTI2018-099655-B-I00, in part by Comunidad de Madrid under Grants
IND2017/TIC-7618, IND2018/TIC-9649, Y2018/TCS-4705 and B2017/BMD-3740 AGES-CM
2CM, in part by BBVA Foundation under 2018 Grant for Deep-DARWiN project, in
part by ISCIII under Grant PI16/01852 and in part by AFSP under Grant
LSRG-1-005-16. We gratefully acknowledge the support of NVIDIA Corporation
with the donation of the Titan X Pascal GPU used for this research.
(Corresponding author: Ignacio Peis).

Ignacio Peis, Pablo M. Olmos, Constanza Vera-Varela, María Luisa Barrigón,
Philippe Courtet, Enrique Baca-García and Antonio Artés-Rodríguez I. Peis, P.
M. Olmos and A. Artés-Rodríguez are with the Dept. of Signal Theory and
Communications, UC3M, Leganés, Madrid, 28911 Spain and also with the Health
Research Institute Gregorio Marañón, 28007 Madrid e-mails<EMAIL_ADDRESS><EMAIL_ADDRESS>antonio@tsc.uc3m.es).C. Vera-Varela, M. L. Barrigón and E.
Baca are with the Dept. of Psychiatry, IIS-Jimenez Díaz Foundation, Madrid,
28040 Spain e-mails<EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS>Courtet is with the Dept. of Psychiatric Emergency and
Acute Care, Lapeyronie Hospital, University of Montpellier, Montpellier, 34295
France e-mail<EMAIL_ADDRESS>L. Barrigón and E. Baca are
also with the Dept. of Psychiatry, Universidad Autónoma de Madrid, Madrid,
28049 Spain.A. Artés-Rodríguez and E. Baca-García are also with the CIBERSAM,
ISCIII, Madrid, 28029 Spain.E. Baca-García is also with Universidad Católica
del Maule, Talca, Chile.

###### Abstract

This article presents a novel method for predicting suicidal ideation from
Electronic Health Records (EHR) and Ecological Momentary Assessment (EMA) data
using deep sequential models. Both EHR longitudinal data and EMA question
forms are defined by asynchronous, variable length, randomly-sampled data
sequences. In our method, we model each of them with a Recurrent Neural
Network (RNN), and both sequences are aligned by concatenating the hidden
state of each of them using temporal marks. Furthermore, we incorporate
attention schemes to improve performance in long sequences and time-
independent pre-trained schemes to cope with very short sequences. Using a
database of 1023 patients, our experimental results show that the addition of
EMA records boosts the system recall to predict the suicidal ideation
diagnosis from 48.13% obtained exclusively from EHR-based state-of-the-art
methods to 67.78%. Additionally, our method provides interpretability through
the t-SNE representation of the latent space. Further, the most relevant input
features are identified and interpreted medically.

###### Index Terms:

Deep learning, RNN, attention, EMA, suicide.

##  I Introduction

Although research on suicidal behavior has intensified in the last decade and
there is now a greater understanding of the factors at play, suicide remains a
complex public health problem [1] and the global prevalence of suicide
continues to be high, with an annual global age-standardized suicide rate of
11.4 per 100 000 people [2]. Globally, lifetime prevalence rates are
approximately 9.2% for suicidal ideation and 60% of transitions from ideation
to plan and attempt occur within the first year after ideation onset [3].

The vast majority of the studies in suicide are based on self report measures
of indicator variables as hopelessness, depressive symptoms, sleep disorders
or violence [4]. These traditional methods for assessment have the
disadvantage of considering long window assessments (weeks, months, or years),
which might exclude crucial reports about who is at imminent risk for suicidal
behavior [5]. In suicidology, as well as in other areas of psychiatry,
technological advances have enhanced ambulatory monitoring techniques and
facilitate the use of Ecological Momentary Assessment (EMA) in combination
with traditional Electronic Health Records (EHR). EMA allows to measure
variables of interest in real-time and in natural contexts of daily life [6],
eliminating inaccurate reporting of the time elapsed between the prompt and
the completion of the questionnaire, and has demonstrated to be safe in
suicidal thoughts assessment [7]. The employment of EMA using smartphones is
opening up new perspectives in suicidal ideation studies, making possible the
identification of different digital phenotypes of suicidal thoughts [8, 9].

In the context of computational psychiatry, there exist a small group of
studies that apply machine learning techniques to predict suicidal behavior.
In [10], the authors develop hand-crafted patient features from longitudinal
EHR data of mental health and primary care visits from a group of nearly 3
million patients to predict suicide attempts and suicide deaths within 90
days, obtaining 0.851 and 0.861 of ROC (Receiver Operating Characteristic),
respectively. In [11], a similar setup is used to predict suicides and
suicidal behavior from longitudinal EHR within three to four years in advance
(0.45 recall / 0.9 precision). In both cases, suicide prediction is performed
using a long-temporal window (3 months and up to 4 years respectively).
Combined with a population of patients in which each of them is followed for
years (in [10], 56% of the patients in the study are treated for more than 5
years, and in [11] the average period is 5.27 years), these results in a well-
condition prediction problem, easily addressed by LASSO (Least Absolute
Shrinkage and Selection Operator) logistic regression [10] and naive Bayesian
classifiers [11]. Clearly, this setup is far from the reality that a clinician
faces, as the longitudinal EHR data for each patient is in most cases less
than one year and a half. In this study, the average periods of the treatments
are 1.3 $\pm$ 0.86 years for EHR and 0.33 $\pm$ 0.44 years for EMA.

In this work, we develop sequential models based on neural-networks [12] for
predicting the identification of suicidal ideation, identified by the
clinician during the presential evaluation at every follow-up visit of the
patient. The difficulty arises in modeling two sources of asynchronous,
length-variable and randomly-sampled sequences, i.e. EHR and EMA. The random
nature of these sequences w.r.t. sampling and length arises in the fact that
they are obtained from patients appointments with the clinician, which are in
general not periodic. Also, the patient’s pattern in accessing EMA mobile
application is in general completely random.

We develop Recurrent Neural Networks (RNN) models able to combine patient’s
EHR longitudinal data along with the EMA question forms. By combining pre-
trained feed-forward Neural Networks (NNs) and attention models [13], our
model learns from the heterogeneous sequences in a fair basis, preventing that
the model training overfits on long temporal sequences for which more
information is available. Using a test patient population, we show our method
is able to reach a 0.67 recall performance, while, when using only EHR
longitudinal data, it downs to 0.46, which is comparable to EHR-based state-
of-the-art methods [10]-[11] that do not explicitly model temporal
dependencies. This demonstrates the great impact that EMA can have in patient
assessment. Further, our model only effectively captures the EMA influence on
the diagnosis when is combined with neural attention models [14]. Therefore,
from both the technical method proposed and the application sides, we believe
this work represents a step forward in the use of machine-learning methods in
suicidology.

The paper is organized as follows. Firstly, Section II introduces the problem
we are facing with (II-A) and the solution we propose (II-B). In Section III
we expose an in-depth presentation of the model architecture. Section IV shows
numerical results and finally Section V contains some concluding remarks.

##  II Problem Statement

In this work we aim at predicting the identification of suicidal ideation by a
mental health clinician (psychiatrist or clinical psychology) using multiple
information sources that contain relevant information about the evolution of a
patient along his/her treatment. On one hand, we use EHR filled up by the
clinicians at every follow-up visit with the patient, which we call _clinician
sequence_. For the EMA data, the so-called _patient sequence_ , we use _MEmind
Wellness Tracker_ 111The web application is available at www.memind.net using
internet‐connected devices with any operating system., a web application for
mobile devices developed to merge different data sources and provide summaries
of the patient state in clinical practice [6].

###  II-A Case Definition

The EHR data variables are organized in four categories: 1) sociodemographic
information; 2) diagnoses according International Statistical Classification
of Diseases and Related Health Problems 10th revision (ICD-10 [15]), questions
regarding history of aggression and violence, suicidal thoughts and behavior
[16], the Clinical Global Impression (CGI) scale that rates the psychiatrist’s
impression of the severity of psychopathology and the variation since last
visit on a similar seven-point scale [17], and the Global Assessment of
Functioning (GAF) Scale [18]; 3) Brugha’s life threatening events scale [19]
and 4) collecting treatment.

The EMA data variables belong to three categories 1) the How are you today?, a
clinical-based questionnaire; 2) the WHO-5 Well-Being Index (World Health
Organization [20]), and 3) the 12-Item General Health Questionnaire (GHQ-12)
[21].

Note that the EHR view also contains the variable we aim at predicting, i.e.,
whether the clinician considers that patient presents _suicidal ideation_ or
not. This variable is not included in the _clinician sequence_ , as it is the
target to be predicted at every visit.

Figure 1: Clinician and patient sequences associated to three different
patients.

We have collected a database formed by the two different sources (EHR and
EMA). The patients were evaluated in psychiatric routine or psychological
visits at mental health facilities affiliated with the Fundación Jiménez Díaz
Hospital in Madrid, Spain. We considered only those patients that accessed the
EMA interface on at least three occasions. From May 2014 to Sep 2016, a sample
of 1023 patients accomplished this criteria. Participants were mostly women
(64.7%), with a mean age of 43.32 $\pm$ 12.6 years old, an active job status
(51.4%) and concerning civil status, 57.2% were married. Most patients were
diagnosed with anxiety related disorders (F40-F49 ICD-10) (52.6%) and with
mood disorders (F30-F39 ICD-10) (23.3%).

Temporal patterns are extremely diverse. For instance, in Fig. 1 we include an
example of the two sequences associated to three patients. Each point
corresponds to an observation. The prediction model must be able to learn from
variable-length sequences in which the typical sequence length is very short
(less than 10 in both cases), but at the same time a certain subset of
patients are characterized with very large sequences. The maximum length of
the patient sequences is 858, while the minimum is 3. The maximum length of
the clinician sequences is 119, while the minimum is only one. The mean
lengths are 19.66 answers and 7.87 visits, respectively.

###  II-B Method proposed

We build a temporal classifier that aims at predicting suicidal ideation at
every visit between the clinician and the patient. Unlike most of the studies
in the related literature, we do take into account the temporal correlation of
the different sequences by using RNNs [12, 22]. A RNN is able to capture
signal correlation by sharing the NN weights along time. Deep Learning
architectures are able to handle the heterogeneity present in the data and
implicitly adapt the variable length of the sequences, as well as capture both
short and long-term temporal correlations. In order to effectively model the
two asynchronous sequences, we propose to use a different RNN for each of
them, where the time index is determined by every follow-up visit in the
clinician sequence, and every time the patient fills up a MEmind questionare
in the patient sequence.

Every event comes with a time stamp that will be used to align both sequences.
As illustrated in Fig. 2, to predict the target variable at every visit, we
concatenate the output of the RNN modeling the clinician sequence with the
most-recently output of the RNN modeling the patient sequence, and the
resulting vector will be used as input to a final classifier built over a
feedforward NN with sigmoidal output activation. In order to limit the number
of parameters and control overfitting, NN layers for dimensionality reduction
are employed. A detail description of the model and its optimization will be
provided in Section III.

As both the clinician RNN and the patient RNN are driven by events and not by
a natural temporal scale, the elapsed time between events will be lost unless
we explicitly encode it as input information to the RNN. Despite RNNs
naturally handle variable-length sequences, the prediction model is trained
over a non-uniform database where more than one-third of the sequences are not
really temporal sequences, as their length is only one. As shown later, a
naive RNN implementation would not result in accurate results. Furthermore,
the use of attention schemes (see Section III-E) is the key to effectively
incorporate the patient sequence. In words, attention schemes learn to
represent variations in the patient state, which is what certainly clinicians
observe, rather than punctual patient MEmind answers.

Figure 2: Model architecture, showing the non-linear transformation of the
inputs and the the merging of the two sequences, adding additional static
information. The figure represents an example of how the model makes a
prediction at the $i$-th patient-clinician visit. Circles refer to
observations, whereas rectangles refer to states.

##  III Design methodology

In this section we present the proposed model, exposed in Figure 2. Firstly,
we describe our notation. Secondly, every part of the architecture is
analyzed.

###  III-A Notation

We deal with EHR (clinician) and EMA (patient) data sources from a group of
$N$ patients. If we use $k_{c}$ and $k_{p}$ variables from the clinician and
patient data, respectively, the $n$-th patient can be represented by two
sequences of $T_{n}^{c}$ and $T_{n}^{p}$ tuples:
$(t_{n,i}^{c},\textbf{x}_{n,i}^{c})\in\mathbb{R}\times\mathbb{R}^{k_{c}},i=1,...,T_{n}^{c}$
and
$(t_{n,j}^{p},\textbf{x}_{n,j}^{p})\in\mathbb{R}\times\mathbb{R}^{k_{p}},j=1,...,T_{n}^{p}$,
where $t_{n,i}^{c}$ and $t_{n,j}^{p}$ denote the time of the $i$-th visit and
$j$-th answer of the $n$-th patient. Each $\textbf{x}_{n,i}^{c}$ and
$\textbf{x}_{n,j}^{p}$ are multivariate observations. We also include the
static (non sequential) information of the patient as $\textbf{x}_{n}^{s}$.
The goal is to predict the binary label $y_{n,i}$ at each clinician time step
$t_{n,i}^{c}$. In the following, we do not include the index of the patient
$n$ for simplicity.

###  III-B Sequence Modeling and Dimensionality Reduction

Provided that we have a fairly small database, and with the aim at avoiding
overfitting, we use basic RNNs instead of more complex architectures like Long
Short Term Memory (LSTM [23]) or Gated Recurrent Unit (GRUs [24]). The RNN is
defined by the following recursion [12]:

$\textbf{h}_{i}=\sigma\left(W^{hx}\textbf{x}_{i}+W^{hh}\textbf{h}_{i-1}+\textbf{b}_{h}\right)$
(1)

where $\sigma$ is a non linear activation (typically hyperbolic tangent or
Rectified Linear Unit), $\textbf{x}_{i}$ the input (in our case, a
multivariate observation with index $i$), $\textbf{h}_{i}$ the so-called RNN
state (with index $i$), $W^{hx}$, $W^{hh}$ projection matrices and
$\textbf{b}_{h}$ the bias term. The RNN output is computed as follows:

$\hat{\textbf{o}}_{i}=\text{softmax}(W^{yh}\textbf{h}_{i}+\textbf{b}_{y})$ (2)

Further, a dimensionality reduction stage is performed over the input data to
both the patient and clinician sequence. To this end, we compared two
different solutions. On one hand, we can use a simple linear layer consisting
on an embedding matrix [25] with dimensions $N\times I$ ($N$ refers to the
dimension of the data, $I$ to the number of inputs for the RNN). This matrix
is trained along with the rest of parameters of the structure in Fig. 2. On
the other hand, we propose to use a nonlinear dimensionality reduction driven
by a dense layer with hyperbolic tangent output activation. This helps us to
better accommodate very short sequences, as this layer is pre-trained
separately over a prediction model in which we neglect temporal correlations
(no RNNs are used) and suicidal ideation $y_{i}$ is predicted directly from
the output of these layers, separately for the clinician (using
$\textbf{x}_{i}^{c}$) and patient (using $\textbf{x}_{j}^{p}$) sources. Due to
the use of these pre-trained dense layers in Fig. 2, the inputs of the RNN are
already discriminative, and the RNN has the specific role of capturing
temporal correlations.

###  III-C Sequence merging and Cost Function

Both the clinician and patient RNNs are aligned in time upon prediction. To
this end, we concatenate the output vector $\textbf{o}_{i}^{c}$ of the
clinician RNN, the selected output $\textbf{o}_{j}^{p}$ (index $j$ of the most
recent answer) of the patient RNN and the static information $\textbf{x}^{s}$,
to form a new vector that is going to feed another feed forward NN that
performs the detection (see Fig. 2). The static information is a 2-dimesional
vector containing a binary feature with the sex of the patient and a scalar
feature with his/her age. The classifier consists on a feed-forward NN that
applies a sigmoid activation at the last layer to obtain predictions between 0
and 1. As described in Section IV, after cross validating several structures,
a simple single dense layer is selected to implement this classifier.

We are using a method similar to the target replication technique exposed in
[26], in which the cost function is given by the sum of the cross entropy
function between the probabilities $\hat{y}_{n,i}$ and the target $y_{n,i}$ at
each $t_{n,i}^{c}$, following:

$\mathcal{L}(\boldsymbol{\Theta})=-\sum_{n=1}^{N}\frac{1}{T_{n}^{c}}\sum_{i=1}^{T_{n}^{c}}\left(y_{n,i}\log\hat{y}_{n,i}+(1-y_{n,i})\log(1-\hat{y}_{n,i})\right)$
(3)

###  III-D Learning the initial state of the RNNs

To reduce dependency w.r.t. initialization, we incorporate a feed forward NN
that, using the initial vector of the clinician sequence and the patient
static information, learns the initial state of the clinician RNN (Figure 3),
and thus we expect an improved adaptation of the model for short sequences.

Figure 3: Schematic of the initial state learning method.

###  III-E Attention schemes

The use of attention schemes [27, 28] within sequence RNN-based models have
become a recent revolution in fields such as Natural Language Processing [13],
or computer vision [29]. In a nutshell, attention models allow neural sequence
models to find discriminative representations in the evolution of the state of
the RNN, applying different weights in the more recent history that can
facilitate learning mid and long-range dependencies. This idea has been
adapted to EHR data in architectures like RETAIN [25], where authors use a
two-level neural attention mechanism with the aim at increasing the model
interpretability. We propose a windowed version of attention scheme based on
[13], where we can tune the size of the past memory. This architecture will
lead us to focus on the most recent answers (see Figure 4) in the patient
sequence, transforming the outputs of the RNN into a new attention-based
space. We take the previous $L$ output vectors of the patient’s RNN (for
sequences shorter than $L$, we apply zero padding) as memory
$\textbf{Y}_{i}=\left[\textbf{o}_{i-L}\,...\,\textbf{o}_{i-1}\in\mathbb{R}^{k\times
L}\right]$, where $k$ is the dimension of the outputs. In figure 4 a schematic
of the proposed attention schemes is included. We define
$\textbf{1}\in\mathbb{R}^{L}$ as a vector of ones. Then, as in [13], we use:

$\displaystyle\textbf{M}_{i}=\tanh\left(\textbf{W}^{Y}\textbf{Y}_{i}+\left(\textbf{W}^{o}\textbf{o}_{i}\right)\textbf{1}^{T}\right)$
$\displaystyle\in\mathbb{R}^{k\times L},$ (4)
$\displaystyle\boldsymbol{\alpha}_{i}=\text{softmax}(\textbf{w}^{T}\textbf{M}_{i})$
$\displaystyle\in\mathbb{R}^{1\times L},$ (5)
$\displaystyle\textbf{r}_{i}=\textbf{Y}_{i}\boldsymbol{\alpha}^{T}$
$\displaystyle\in\mathbb{R}^{k}$ (6)

where $\textbf{W}^{Y}$, $\textbf{W}^{o}\in\mathbb{R}^{k\times k}$ and w are
trainable parameters. Vector $\textbf{r}_{i}$ is the attention-weighted
representation of the previous outputs. Thus, the transformed output at time
instant with index $i$ is given by:

$\displaystyle\textbf{o}_{i}^{*}=\tanh\left(\textbf{W}^{r}\textbf{r}_{i}+\textbf{W}^{x}\textbf{o}_{i}\right)$
$\displaystyle\in\mathbb{R}^{k}$ (7)

where $\textbf{W}^{r}$, $\textbf{W}^{x}\in\mathbb{R}^{k\times k}$ are also
trainable projection matrices.

The experimental results included in the next section show that attention
mechanisms provide a dramatic performance improvement even when dealing with
moderately long (length above 4) sequences, which in turn compose the major
part of our dataset.

Figure 4: Proposed model architecture, including an attention scheme on the
patient sequence. For this graphical example, we use $L=2$.

##  IV Numerical Experiments

The results are obtained using a total of 768 patients for training and 255
for test. Given the nature of the problem to solve, we are more interested in
evaluating the recall. K-fold validation is used to validate structural
parameters, regularization parameters, and the number of epochs (early
stopping). To reduce inter-fold variance, training set in each fold is
characterized by roughly a uniform distribution of EHR sequences lengths.
I.e., each fold contains the same proportion of sequences of lengths 1, 2, 3,
and 4 or more. In Table I all the validated parameters for the proposed
architecture are included. Following early stopping technique, 100 epochs are
chosen for preventing overfitting. Model regularization is achieved by
combining $L_{2}$ regularization of all the weights in the structure and
dropout at the input of the final classifier [30].

Parameter | Value
---|---
RNN hidden units | [clinician: 10, patient: 5]
Input NN hidden units | [layer_1: 10, layer_2: 20]
Output NN hidden units | [layer_1: 10]
Initial state NN hidden units | [layer_1: 10, layer_2: 10]
Dropout probability | 0.6
Epochs | 100
Bath size | 20
Train size | 0.75
Learning rate | 0.005
L2 regularizer factor | 0.01

TABLE I: Parameter configuration for the model architecture.

###  IV-A Performance discussion

| Train | Test
---|---|---
Length | All | $1$ | $2$ | $3$ | $\geq 4$ | All
NN | 1.14 $\pm$ 0.37 | 0.00 $\pm$ 0.00 | 0.00 $\pm$ 0.00 | 0.00 $\pm$ 0.00 | 0.70 $\pm$ 1.22 | 0.45 $\pm$ 0.43
LogReg | 65.51 $\pm$ 6.06 | 37.50 $\pm$ 41.46 | 35.41 $\pm$ 20.73 | 43.75 $\pm$ 44.63 | 51.29 $\pm$ 5.12 | 51.24 $\pm$ 5.23
HMM | 59.92 $\pm$ 6.70 | - | 2.08 $\pm$ 3.62 | 12.50 $\pm$ 21.65 | 51.00 $\pm$ 10.48 | 52.38 $\pm$ 14.62
Linear inputs | 76.43 $\pm$ 5.98 | 25.00 $\pm$ 43.30 | 16.67 $\pm$ 20.41 | 31.25 $\pm$ 44.63 | 58.64 $\pm$ 8.14 | 55.22 $\pm$ 5.66
Non-linear inputs | 77.30 $\pm$ 3.69 | 25.00 $\pm$ 43.30 | 25.00 $\pm$ 21.65 | 43.75 $\pm$ 44.63 | 57.78 $\pm$ 4.93 | 56.14 $\pm$ 4.73
Attention, L=1 | 85.16 $\pm$ 9.92 | 25.00 $\pm$ 43.30 | 25.00 $\pm$ 25.00 | 43.75 $\pm$ 44.63 | 68.97 $\pm$ 4.08 | 67.68 $\pm$ 3.50
Attention, L=2 | 80.70 $\pm$ 3.55 | 25.00 $\pm$ 43.30 | 25.00 $\pm$ 25.00 | 43.75 $\pm$ 44.63 | 64.45 $\pm$ 7.16 | 63.59 $\pm$ 7.10
Attention, L=3 | 81.63 $\pm$ 5.02 | 25.00 $\pm$ 43.30 | 25.00 $\pm$ 25.00 | 43.75 $\pm$ 44.63 | 66.69 $\pm$ 3.60 | 65.57 $\pm$ 3.12
Attention, L=4 | 80.76 $\pm$ 7.47 | 25.00 $\pm$ 43.30 | 25.00 $\pm$ 25.00 | 43.75 $\pm$ 44.63 | 72.63 $\pm$ 7.53 | 71.16 $\pm$ 7.04

TABLE II:  Recall provided by different models, including variations of our
proposal. We test the model for every group of EHR sequences (depending on the
length). All the values correspond to the mean and standard deviation of the
metric obtained with 4 folds. Decision threshold is selected for each case by
cross-validation.

Table II shows recall average performance and standard deviation measured
along the different folds of three baseline models and five different
variations of our model. Test performance is stratified, considering clinician
EHR sequences of length 1, 2, 3, and 4 or more. However, we note that the vast
majority of suicidal ideation cases (more specifically, 65.43% of true
positives) in the dataset are contained in long EHR sequences (4 or more).
Namely, achieving a good performance in the last group is crucial.

With NN we refer to a feed-forward neural network of two layers with 10 and 20
hidden units, respectively. LogReg corresponds to a simple logistic regression
with $L_{2}$ norm regularization, which is one of the simplest classifiers.
These two cases are trained using the stacked vectors of visits
$\textbf{x}_{n,i}^{c}$ and the most recent answers $\textbf{x}_{n,j}^{p}$,
i.e., sequential nature is ignored. HMM is a first order Hidden Markov Model
with both Gaussian and labels emission and 5 hidden states.

With the fourth and fifth rows we can compare the effect in our method of
using a linear transformation, versus a pre-trained dense NN layer to perform
dimensionality reduction at the RNN input (see Section III-B). It can be
observed that the use of a pre-trained dense NN layer, which already provides
discriminative features to the RNN, provides both better train and test
performance, particularly for short sequences. The last four rows correspond
to the model versions that use attention-based patient RNN outputs, with
attention windows varying from $L=1$ to $L=4$. The version with $L=1$ is
chosen as the most suitable, considering a trade-off between recall
performance, computational cost, and results’ standard deviation. Note that
for $L=1$ the patient RNN contains features related to the variation between
two consecutive RNN states. Observe the remarkable performance improvement in
this latter case, particularly for long clinician sequences, which correspond
to patients that have been under psychiatric treatment for a sufficiently long
time.

For short EHR sequences, all models provide poor recall and large inter-fold
standard deviation. This is certainly related to the lack of data of this type
in our dataset. Note however that the use of non-linear RNN inputs in our
model helps to cope with EHR sequences with length equal to two. Even without
attention, our proposal clearly improves the baselines (NN, LogReg, and HMM)
for long (length 3 or more) EHR sequences. Both Logistic Regression and HMM
are comparable (they both achieve roughly a recall of 51% for long sequences),
while NN certainly is not able to predict the suicidal ideation cases. Our
model with attention and $L=1$ at the patient sequence achieves a recall close
to 67% with the smallest variance across folds. For this model, other metrics
evaluations are included in table III for the chosen version.

| Train | Test
---|---|---
Length | All | $1$ | $2$ | $3$ | $\geq 4$ | All
Accuracy | 95.98 $\pm$ 0.60 | 92.65 $\pm$ 9.15 | 91.26 $\pm$ 5.91 | 91.00 $\pm$ 9.38 | 88.76 $\pm$ 5.11 | 88.96 $\pm$ 5.04
Precision | 87.08 $\pm$ 6.63 | 25.00 $\pm$ 43.30 | 33.18 $\pm$ 34.82 | 19.33 $\pm$ 21.05 | 58.88 $\pm$ 5.11 | 58.81 $\pm$ 20.06
Recall | 85.16 $\pm$ 9.92 | 25.00 $\pm$ 43.30 | 25.00 $\pm$ 25.00 | 43.75 $\pm$ 44.63 | 68.97 $\pm$ 4.08 | 67.68 $\pm$ 3.50
AUC | 94.02 $\pm$ 4.22 | 58.65 $\pm$ 38.40 | 67.33 $\pm$ 34.82 | 59.07 $\pm$ 37.48 | 83.85 $\pm$ 21.05 | 83.29 $\pm$ 5.13
F1-score | 85.34 $\pm$ 4.02 | 25.00 $\pm$ 43.30 | 28.43 $\pm$ 28.62 | 26.11 $\pm$ 26.68 | 61.83 $\pm$ 11.81 | 61.33 $\pm$ 4.03

TABLE III: Performance metrics of the model proposed with $L=1$ attention.

The results in Table II demonstrate that EMA is a significant support for
automatically detect suicidal ideation when is used along with the clinician
information. Results indeed demonstrate that inspecting the contrast between
two consecutive questionnaires rather than basing the prediction exclusively
on the latest sets of responses improves the performance. This information is
exactly what the attention model is attempting to capture. A similar
conclusion is drawn when we compare the performance of the model when either
clinician RNN, or the patient RNN are not used (all parameters are set to 0
and the RNN is not trained so it does not influence the classifier). Fig. 5
shows the recall achieved by these models as a function of the classifier
threshold $\alpha$. Not surprisingly, when the clinician RNN is removed and
only the EMA information is used, the model is no longer valid to predict
suicidal ideation diagnose, as the EMA sequence is not aligned in time and
completely independent of the suicidal ideation identification sequence.
However, when both sources of information are used and attention is used to
process the patient RNN prediction recalls improves in close to a 20% w.r.t.
to the model when the patient RNN is not used or only the latest Memind
questionnaire is considered.

Figure 5: Test recall results evaluating the merging effect of the sequences,
depending on the threshold for detection. The lines represent the results
obtained using only the patient (blue) or the clinician (orange) sources,
using both (green) and using both with attention schemes (red).

###  IV-B Most relevant features

After model has been validated in terms of test performance, we can study what
input features in both the clinician sequence and patient sequence impact the
most in the classifier’s prediction. To this end, we analyze the profile of
the weights in the first layer of the NN that performs the input
dimensionality reduction. Thus, we have a set of 93 values associated to each
of the clinician input features and 24 for the patient features. Tables IV and
V shows the top ten for the clinician and patient sources, respectively. We
also indicate the associated normalized weight (absolute value of the weight
divided by the maximum). Observed that, while in the clinician sequence we
have four input features that are clearly dominant, a similar effect is not
present in the patient sequence. This is indeed expected, since we have
observed that the patient sequence really matters when two consecutive
questionnaires are jointly processed using the attention mechanism. It is
important to highlight that two questions out of the total of five directly
related to suicidal risk were not found within the dominant input features of
the clinician sequence, these are suicidal plans and suicide attempt since the
last visit. This finding is consistent with previous studies that criticize
the traditional temporal sequence of stages in the suicidal process [31].
Furthermore, it is interesting to remark how diagnoses less commonly
associated with suicidal ideation such as F00-F09 diagnoses are found to be of
greater importance compared to those most frequently studied (F40-F49 and
F30-F39 ICD-10 codes). Although not shown in Table III, it is worth to mention
that the 12-th question among the largest weights corresponds to the
unemployed item in Brugha’s life threatening events scale, which is a common
clinical finding. In EMA sequence, questions related to despair, aggression,
confidence or rest have largest weights than more specific questions like
“Have you ever felt that you had no desire to live?”, which is in number 13 of
the 24 input features. This is of clinical interest, as indirect and non-
intrusive questions can make clinicians be alert about suicidal thoughts in
patients.

Clinician
---
Index | Score | Feature
1 | 1.00 | Have you ever thought about harming yourself?
2 | 0.81 |  | Have you ever thought you would be better
---
off dead or wish you were dead?
3 | 0.65 | Have you ever attempted suicide?
4 | 0.51 | Unemployed with benefit.
5 | 0.48 | Living with siblings.
6 | 0.45 | Shared residence.
7 | 0.44 |  | F00-F09 Mental disorders due to known
---
physiological conditions.
8 | 0.41 | Electroconvulsive therapy.
9 | 0.40 | Living with offspring.
10 | 0.40 | Temporary disability.

TABLE IV: Top ten more relevant features from the clinician (EHR) source.
Patient
---
Index | Score | Feature
1 | 1.000 |  | Have you recently felt you could no overcome your
---
difficulties? (GHQ-12, item 6)
2 | 0.98 | Anger, arguments or fights.
3 | 0.92 |  | Have you recently been thinking of yourself as a
---
worthless person? (GHQ-12, item 11)
4 | 0.91 |  | I have felt active and vigorous (WHO-5,
---
item 3)
5 | 0.91 |  | Have you recently felt capable of making decisions
---
about things (GHQ-12, item 4)
6 | 0.90 |  | Have you recently been losing confidence in
---
yourself? (GHQ-12, item 10)
7 | 0.90 |  | Have you recently been feeling reasonably happy,
---
all things considered? (GHQ-12, item 12)
8 | 0.89 |  | I woke up feeling fresh and rested (WHO-5,
---
item 4)
9 | 0.89 |  | Have you recently felt constantly under strain?
---
(GHQ-12, item 5)
10 | 0.88 |  | Have you recently been able to enjoy your normal
---
day to day activities? (GHQ-12, item 7)

TABLE V: Top ten more relevant features from the patient (EMA) source.

###  IV-C Latent Space Visualization

One of the main advantages of deep learning methods is the ability of
implicitly finding latent representation of the data. For that purpose we are
employe t-SNE [32] (t-Distributed Stochastic Neighbor Embedding), a powerful
technique for visualizing high-dimensional data using a transformation that
reveals structure at many different scales. We are applying this technique on
the vectors that feed the final classifier to evaluate the capacity of
discrimination of the RNN networks. For all the experiments the t-SNE
perplexity parameter is fixed to 100 and the learning rate to 10. In Figure 6
we include two examples of t-SNE visualizations of the outputs of the RNNs. At
the left we include the observed suicidal ideation corresponding to each point
(i.e. the label we aim at predicting), while at the right we plot the
probability of suicidal ideation that the model obtains. Indeed, we include
trajectories in the map from three different patients. It seems clear than the
t-SNE transformation reveals clouds of points with groups of patients that
present suicidal ideation, and the state ‘travels’ between them depending on
the evolution of the suicidal ideation of the patient. As a by product,
observe that in the right map all positives are concentrated and do not
overfit to the isolated positive targets that can be seen in the left map.
This is a sign that the model is not overfitting the training data.

(a) Real suicidal ideation (b) Suicidal ideation probability

Figure 6: Dimensionality reduction of the latent space (the vector that
concatenates $\textbf{o}_{i}^{c}$, $\textbf{o}_{j}^{p}$ and $\textbf{x}^{s}$
in Figures 2 and 4), using t-SNE. Each point corresponds to a visit. Left:
original suicidal ideation. Right: the color map is the probability of
suicidal ideation (obtained logits of the final feed forward network). Indeed,
trajectories in the map for three different patients are included.

##  V Discussion

A model for asynchronous, length-variable, randomly sampled sequences with
heterogeneous data has been developed, using RNNs for the sequential inputs
with a first non linear transformation, given by pre-trained feed-forward NNs,
to cope with the heterogeneity of the data and perform a dimensionality
reduction. The analysis of the input layer of the feed-forward networks has
allowed us to obtain a measure of the features relevance, and perform clinical
interpretations. Furthermore, attention mechanisms are incorporated to the
_patient sequence_ for addressing the variability in sequence lengths. We have
applied a merging technique consisting on temporally aligning the outputs of
each RNN according to time labels. This merged vectors feed another feed
forward NN that predicts the suicidal ideation at each visit to the clinician.
In our case, the suicidal ideation sequence imposed the date and we extracted
from each RNN the output at last time instant before that date. We have built
a cost function that takes the sum of the cost in all visits of a patient. We
have demonstrated that this merging procedure improves the performance of the
model, rather than if each sequence was used separately.

Our model predicts suicidal ideation interpreted by a clinician, based on the
past mental health visits and the history of answers of a patient collected
through a web-based EMA system, that provides a promising complementary tool
for the clinician. The visualization step has allowed us to discover the
trajectories of the patient within the latent space, observing that there
exist some clouds with positive and negatives samples, and there exists a
coherence between the sequential target of the patient and the trajectory on
the t-SNE map.

This is, to our knowledge, a novel study with two main innovations: 1 the
combination of Electronic Health Records and Ecological Momentary Assessment
data for predicting suicidal ideation, and 2 the use of Deep Learning
sequential architectures (RNNs) to discover temporal correlations in
psychiatric patients’ evolution. We have demonstrated that our method
outperforms other existing methods in the literature. Further, we have exposed
an interesting feature of our model: the interpretability of the latent space
via dimensionality reduction techniques as t-SNE.

Future work could be centered on including more sequences associated to the
patient that can be merged and aligned, for instance, the complete Electronic
Health Record of the patient (including primary care visits). Another branch
of the research will be the application of more complex attention schemes on
every sequence, with different temporal windows. At last, we will research in
such an interesting application of this work consisting on generative models
for this medical data.

## References

  * [1] G. Turecki and D. A. Brent, “Suicide and suicidal behaviour,” _The Lancet_ , vol. 387, no. 10024, pp. 1227–1239, 2016.
  * [2] Suicide prevention. [Online]. Available: http://www.who.int/mental_health/suicide-prevention/en/
  * [3] M. K. Nock, G. Borges, E. J. Bromet, J. Alonso, M. Angermeyer, A. Beautrais, R. Bruffaerts, W. T. Chiu, G. De Girolamo, S. Gluzman _et al._ , “Cross-national prevalence and risk factors for suicidal ideation, plans and attempts,” _The British Journal of Psychiatry_ , vol. 192, no. 2, pp. 98–105, 2008.
  * [4] C. L. Davidson, M. D. Anestis, and P. M. Gutierrez, “Ecological momentary assessment is a neglected methodology in suicidology,” _Archives of suicide research_ , vol. 21, no. 1, pp. 1–11, 2017.
  * [5] E. Czyz, C. King, and I. Nahum-Shani, “Ecological assessment of daily suicidal thoughts and attempts among suicidal teens after psychiatric hospitalization: Lessons about feasibility and acceptability,” _Psychiatry Research_ , 2018\.
  * [6] M. L. Barrigón, S. Berrouiguet, J. J. Carballo, C. Bonal-Giménez, P. Fernández-Navarro, B. Pfang, D. Delgado-Gómez, P. Courtet, F. Aroca, J. Lopez-Castroman _et al._ , “User profiles of an electronic mental health tool for ecological momentary assessment: Memind,” _International journal of methods in psychiatric research_ , vol. 26, no. 1, p. e1554, 2017.
  * [7] M. Husky, E. Olié, S. Guillaume, C. Genty, J. Swendsen, and P. Courtet, “Feasibility and validity of ecological momentary assessment in the investigation of suicide risk,” _Psychiatry research_ , vol. 220, no. 1-2, pp. 564–570, 2014.
  * [8] E. M. Kleiman, B. J. Turner, S. Fedor, E. E. Beale, R. W. Picard, J. C. Huffman, and M. K. Nock, “Digital phenotyping of suicidal thoughts,” _Depression and anxiety_ , 2018.
  * [9] J. Bernanke, B. Stanley, and M. Oquendo, “Toward fine-grained phenotyping of suicidal behavior: the role of suicidal subtypes,” _Molecular psychiatry_ , vol. 22, no. 8, p. 1080, 2017.
  * [10] G. E. Simon, E. Johnson, J. M. Lawrence, R. C. Rossom, B. Ahmedani, F. L. Lynch, A. Beck, B. Waitzfelder, R. Ziebell, R. B. Penfold _et al._ , “Predicting suicide attempts and suicide deaths following outpatient visits using electronic health records,” _American Journal of Psychiatry_ , pp. appi–ajp, 2018.
  * [11] Y. Barak-Corren, V. M. Castro, S. Javitt, A. G. Hoffnagle, Y. Dai, R. H. Perlis, M. K. Nock, J. W. Smoller, and B. Y. Reis, “Predicting suicidal behavior from longitudinal electronic health records,” _American journal of psychiatry_ , vol. 174, no. 2, pp. 154–162, 2016.
  * [12] Z. C. Lipton, J. Berkowitz, and C. Elkan, “A critical review of recurrent neural networks for sequence learning,” _arXiv preprint arXiv:1506.00019_ , 2015.
  * [13] M. Daniluk, T. Rocktäschel, J. Welbl, and S. Riedel, “Frustratingly short attention spans in neural language modeling,” _arXiv preprint arXiv:1702.04521_ , 2017.
  * [14] D. Bahdanau, K. Cho, and Y. Bengio, “Neural machine translation by jointly learning to align and translate,” _arXiv preprint arXiv:1409.0473_ , 2014\.
  * [15] O. mondiale de la santé, W. H. Organization, W. Staff, and WHO, _The ICD-10 classification of mental and behavioural disorders: clinical descriptions and diagnostic guidelines_. World Health Organization, 1992, vol. 1.
  * [16] T. Hergueta, R. Baker, and G. C. Dunbar, “The mini-international neuropsychiatric interview (mini): the development and validation of a structured diagnostic psychiatric interview for dsm-ivand icd-10,” _J clin psychiatry_ , vol. 59, no. Suppl 20, p. 2233, 1998.
  * [17] W. Guy, “Early clinical drug evaluation (ecdeu) assessment manual,” _Rockville: National Institute of Mental Health_ , 1976.
  * [18] R. L. Spitzer, M. Gibbon, J. Williams, and J. Endicott, “Global assessment of functioning (gaf) scale,” _Outcomes assessment in clinical practice_ , pp. 76–78, 1996.
  * [19] T. Brugha, P. Bebbington, C. Tennant, and J. Hurry, “The list of threatening experiences: a subset of 12 life event categories with considerable long-term contextual threat,” _Psychological medicine_ , vol. 15, no. 1, pp. 189–194, 1985.
  * [20] The WHO-5 website. [Online]. Available: https://www.psykiatri-regionh.dk/who-5/Pages/default.aspx
  * [21] M. D. P. Sánchez-López and V. Dresch, “The 12-item general health questionnaire (ghq-12): reliability, external validity and factor structure in the spanish population,” _Psicothema_ , vol. 20, no. 4, 2008.
  * [22] A. Graves, _Supervised Sequence Labelling with Recurrent Neural Networks_. Springer, 2012.
  * [23] S. Hochreiter and J. Schmidhuber, “Long short-term memory,” _Neural computation_ , vol. 9, no. 8, pp. 1735–1780, 1997.
  * [24] K. Cho, B. Van Merriënboer, C. Gulcehre, D. Bahdanau, F. Bougares, H. Schwenk, and Y. Bengio, “Learning phrase representations using RNN encoder-decoder for statistical machine translation,” _arXiv preprint arXiv:1406.1078_ , 2014.
  * [25] E. Choi, M. T. Bahadori, J. Sun, J. Kulas, A. Schuetz, and W. Stewart, “Retain: An interpretable predictive model for healthcare using reverse time attention mechanism,” in _Advances in Neural Information Processing Systems_ , 2016, pp. 3504–3512.
  * [26] Z. C. Lipton, D. C. Kale, C. Elkan, and R. Wetzel, “Learning to diagnose with lstm recurrent neural networks,” _arXiv preprint arXiv:1511.03677_ , 2015\.
  * [27] Y. Deng, Y. Kim, J. Chiu, D. Guo, and A. M. Rush, “Latent alignment and variational attention,” _arXiv preprint arXiv:1807.03756_ , 2018.
  * [28] Y. Kim, C. Denton, L. Hoang, and A. M. Rush, “Structured attention networks,” _arXiv preprint arXiv:1702.00887_ , 2017.
  * [29] V. Mnih, N. Heess, A. Graves _et al._ , “Recurrent models of visual attention,” in _Advances in neural information processing systems_ , 2014, pp. 2204–2212.
  * [30] N. Srivastava, G. Hinton, A. Krizhevsky, I. Sutskever, and R. Salakhutdinov, “Dropout: a simple way to prevent neural networks from overfitting,” _The Journal of Machine Learning Research_ , vol. 15, no. 1, pp. 1929–1958, 2014.
  * [31] E. Baca-Garcia, M. M. Perez-Rodriguez, M. A. Oquendo, K. M. Keyes, D. S. Hasin, B. F. Grant, and C. Blanco, “Estimating risk for suicide attempt: Are we asking the right questions?: Passive suicidal ideation as a marker for suicidal behavior,” _Journal of affective disorders_ , vol. 134, no. 1-3, pp. 327–332, 2011.
  * [32] L. v. d. Maaten and G. Hinton, “Visualizing data using t-SNE,” _Journal of machine learning research_ , vol. 9, no. Nov, pp. 2579–2605, 2008\.