# Formal differential variables and an abstract chain rule

Samuel Alexander

###### Abstract.

One shortcoming of the chain rule is that it does not iterate: it gives the
derivative of $f(g(x))$, but not (directly) the second or higher-order
derivatives. We present iterated differentials and a version of the
multivariable chain rule which iterates to any desired level of derivative. We
first present this material informally, and later discuss how to make it
rigorous (a discussion which touches on formal foundations of calculus). We
also suggest a finite calculus chain rule (contrary to Graham, Knuth and
Patashnik’s claim that “there’s no corresponding chain rule of finite
calculus”).

##  1\. Introduction

Consider the following statement, uncontroversial in an elementary calculus
context ($*$): “For all variables $u$ and $v$, $d(uv)=v\,du+u\,dv$.” In his
popular calculus textbook [11], Stewart says:

> …the differential $dx$ is an independent variable…

So if $*$ really does hold for _all_ variables $u$ and $v$, and if $x$ is a
variable, and if (as Stewart says) $dx$ is also a variable, then, by letting
$u=x$ and $v=dx$, we get $d(x\,dx)=dx\,dx+x\,ddx$. We do not know whether
Stewart intended us to make such an unfamiliar-looking conclusion from his
innocent-looking statement, but let’s continue along these lines and see where
it leads us. We will formalize this kind of computation using machinery from
first-order logic, and show that it leads to an elegant higher-order
multivariable chain rule.

A weakness of the familiar chain rule is that it does not iterate: it tells us
how to find the _first_ derivative of $f(g(x))$, but it does not tell us how
to find second- or higher-order derivatives of the same (at least not
directly). Our abstract chain rule will iterate: the exact same rule which
tells us $df(g(x))$ will also tell us $d^{k}f(g(x))$ for any integer $k>1$.

Our $d$ operator has some similarities with the $\Delta$ operator of Huang et
al [7]. Our work improves on theirs in that we explicitly distinguish
differential variables from others, so that the operator we develop better
reveals the connection to differentials. For example, in Huang et al, one has
$\Delta_{1}e^{x_{0}}=e^{x_{0}}x_{1}$ and
$\Delta_{2}e^{x_{0}}=e^{x_{0}}(x_{1}^{2}+x_{2})$, which is equivalent to our
$de^{x_{0}}=e^{x_{0}}\,dx_{0}$ and
$d^{2}e^{x_{0}}=e^{x_{0}}(dx_{0}\,dx_{0}+ddx_{0})$. Besides better emphasizing
the connection to differentials, the latter version should also be more
familiar, since we already routinely write things like $de^{x}=e^{x}\,dx$ in
elementary calculus classes.

##  2\. Computing iterated partial derivatives: informal examples

In this section, we will informally describe a way to compute iterated partial
derivatives of a multivariable function. We will make the method formal in
subsequent sections.

######  Example 2.1.

Compute the differential $dd\,x^{2}=d(dx^{2})$, treating differential
variables just like ordinary variables.

###### Solution.

The differential $dx^{2}=2x\,dx$ involves two variables: $x$ and $dx$. Thus,
$d(dx^{2})$ will have two terms, one where we differentiate with respect to
$x$ and multiply the result by $dx$, and one where we differentiate with
respect to $dx$ and multiply the result by $ddx$:

$\displaystyle ddx^{2}$ $\displaystyle=d(dx^{2})$ $\displaystyle=d(2x\,dx)$
$\displaystyle=\frac{\partial(2x\,dx)}{\partial
x}\,dx+\frac{\partial(2x\,dx)}{\partial dx}\,ddx$
$\displaystyle=2\,dx\,dx+2x\,ddx.$

Note that when we compute $\frac{\partial(2x\,dx)}{\partial x}$, we treat $dx$
as a variable independent from $x$, so $dx$ can be treated as a constant.
Likewise when we compute $\frac{\partial(2x\,dx)}{\partial dx}$, $x$ is
treated as a constant. ∎

######  Example 2.2.

Compute the differential $dd\,e^{x}$, treating differential variables just
like ordinary variables.

###### Solution.

As in Example 2.1, since $d\,e^{x}=e^{x}\,dx$,

$\displaystyle dd\,e^{x}$ $\displaystyle=d\,(e^{x}\,dx)$
$\displaystyle=\frac{\partial(e^{x}\,dx)}{\partial
x}\,dx+\frac{\partial(e^{x}\,dx)}{\partial dx}\,ddx$
$\displaystyle=e^{x}\,dx\,dx+e^{x}\,ddx.$

∎

######  Example 2.3.

Compute $dd\,f(x)$, treating differential variables just like ordinary
variables.

###### Solution.

Just as above,

$\displaystyle dd\,f(x)$ $\displaystyle=d(f^{\prime}(x)\,dx)$
$\displaystyle=\frac{\partial(f^{\prime}(x)\,dx)}{\partial
x}\,dx+\frac{\partial(f^{\prime}(x)\,dx)}{\partial dx}\,ddx$
$\displaystyle=f^{\prime\prime}(x)\,dx\,dx+f^{\prime}(x)\,ddx.$

∎

In a later section, we will formalize and prove a formal chain rule (Corollary
6.9). For now, we will state it informally:

######  Remark 2.4.

(Abstract Chain Rule, stated informally) Let $T$ and $U$ be expressions and
let $x$ be a non-differential variable. Assume $T$, $U$, and all of their sub-
expressions are everywhere infinitely differentiable. Then

$d(T[x|U])=(dT)[x|U],$

where the operator $[x|U]$ works by simultaneously replacing all occurrences
of $x$ by $U$, all occurrences of $dx$ by $dU$, all occurrences of $d^{2}x$ by
$d^{2}U$, and so on.

The Abstract Chain Rule can be stated in English: “substituting first and then
applying $d$ gives the same result as applying $d$ first and then
substituting, provided that when one substitutes $U$ for $x$, one also
substitutes $dU$ for $dx$ and so on.”

######  Example 2.5.

Compute $(e^{x^{2}})^{\prime\prime}$.

###### Solution.

By Example 2.3, $(e^{x^{2}})^{\prime\prime}$ is the $dx\,dx$-coefficient of
$dd\,e^{x^{2}}$. We compute:

$\displaystyle dd\,e^{x^{2}}$ $\displaystyle=dd\,(e^{x}[x|x^{2}])$
$\displaystyle=(dd\,e^{x})[x|x^{2}]$ (Abstract Chain Rule)
$\displaystyle=(e^{x}\,dx\,dx+e^{x}\,ddx)[x|x^{2}]$ (Example 2.2)
$\displaystyle=e^{x^{2}}\,d(x^{2})\,d(x^{2})+e^{x^{2}}\,dd(x^{2})$
(Substituting)
$\displaystyle=e^{x^{2}}\,(2x\,dx)^{2}+e^{x^{2}}\,(2\,dx\,dx+2x\,ddx)$
(Example 2.1) $\displaystyle=(4x^{2}+2)e^{x^{2}}\,dx\,dx+2xe^{x^{2}}\,ddx.$

The answer is the above $dx\,dx$-coefficient:

$(e^{x^{2}})^{\prime\prime}=(4x^{2}+2)e^{x^{2}}.$

∎

Our Abstract Chain Rule works for multivariable and higher-order derivatives,
too.

######  Example 2.6.

The iterated total derivative

$d^{3}\sin xy=d^{3}(\sin x\,[x|xy])=(d^{3}\sin x)[x|xy]$

encodes:

  * •

$\partial^{3}\sin xy/\partial x^{3}$ as its $dx\,dx\,dx$-coefficient.

  * •

$\partial^{3}\sin xy/\partial y^{3}$ as its $dy\,dy\,dy$-coefficient.

  * •

$\frac{\partial^{3}\sin xy}{\partial x\partial y\partial
y}=\frac{\partial^{3}\sin xy}{\partial y\partial x\partial
y}=\frac{\partial^{3}\sin xy}{\partial y\partial y\partial x}$ times $3$ as
its $dx\,dy\,dy=dy\,dx\,dy=dy\,dy\,dx$-coefficient (the fact that there are
three ways to write this coefficient is why we write “times 3”).

In Sections 5–6 we will formalize and prove the Abstract Chain Rule. But
first, we will connect these higher-order differentials to a more concrete
higher-order chain rule known as Faà di Bruno’s formula, and also show how the
same ideas lead to a finite calculus chain rule.

##  3\. Faà di Bruno’s formula

Faà di Bruno’s formula, named after the 19th century Italian priest Francesco
Faà di Bruno, is a formula for the higher derivatives of $f(g(x))$. See [8]
and [3] for the history of Faà di Bruno’s formula (see also [10] for related
work in category theory by another ACMS presenter). The formula can be stated
combinatorially:

$f(g(x))^{(n)}=\sum_{\pi\in\Pi_{n}}f^{(|\pi|)}(g(x))\prod_{B\in\pi}g^{(|B|)}(x)$

where $\pi$ ranges over the set $\Pi_{n}$ of all partitions of
$\\{1,\ldots,n\\}$ (so for each such partition $\pi$, $B$ ranges over the
blocks in $\pi$).

The ideas of Section 2 offer an intuitive way to understand the above
formula111Shortly after presenting this argument at ACMS, we realized that the
argument can actually be applied directly, without using iterated
differentials at all, yielding a shockingly short elementary proof of Faà di
Bruno’s formula. Examining the literature, we found that the basic idea is
already known [9] [6], but both published proofs which we found are actually
proofs of more complicated multivariable generalizations of Faà di Bruno’s
formula. For the single-variable special case, the idea (essentially the same
idea which we presented using iterated differentials at ACMS) is so simple
that it can be written with a single sentence [2].. For any partition
$\pi=\\{B_{1},\ldots,B_{k}\\}$ of $\\{1,\ldots,n\\}$, let $I(\pi)$ be the
expression

$I(\pi)=f^{(k)}(x)\,d^{|B_{1}|}x\,d^{|B_{2}|}x\,\cdots\,d^{|B_{k}|}x$

involving iterated differentials as in Section 2. By an inductive argument,
one can check that

$d^{n}f(x)=\sum_{\pi\in\Pi_{n}}I(\pi)$

(for the inductive step, consider the different ways of obtaining a partition
$\pi^{\prime}\in\Pi_{n+1}$ from a partition $\pi\in\Pi_{n}$: one can either
add $\\{n+1\\}$ as a new block, which corresponds to changing $f^{(k)}(x)$ to
$f^{(k+1)}(x)dx$ when using the product rule to calculate $dI(\pi)$; or one
can add $n+1$ to existing block $B_{i}$ of $\pi$, which corresponds to
changing $d^{|B_{i}|}x$ to $d^{|B_{i}|+1}x$ when using the product rule to
calculate $dI(\pi)$).

By similar reasoning as in Examples 2.3 and 2.5, $f(g(x))^{(n)}$ is the
$(dx)^{n}$-coefficient of
$d^{n}f(g(x))=d^{n}(f(x)[x|g(x)])=(d^{n}f(x))[x|g(x)]$. Thus $f(g(x))^{(n)}$
is the $(dx)^{n}$-coefficient of

$\sum_{\pi\in\Pi_{n}}I(\pi)[x|g(x)]=\sum_{\pi\in\Pi_{n}}f^{(|\pi|)}(g(x))\prod_{B\in\pi}d^{|B|}g(x).$

One can check that $d^{|B|}g(x)=g^{(|B|)}(x)d^{|B|}x+o$ where $o$ is a sum of
terms involving higher-order differentials (which can be ignored because they
contribute nothing to the $(dx)^{n}$-coefficient we seek). Faà di Bruno’s
formula follows.

##  4\. Application to finite calculus

The ideas in this paper also lead to a chain rule for the so-called finite
calculus. The finite calculus is described in Section 2.6 of Graham, Knuth and
Patashnik’s _Concrete Mathematics_ [4]. In finite calculus, one defines an
operator $\Delta$ on functions by $\Delta f(x)=f(x+1)-f(x)$. This operator has
many surprising analogies with differentiation, but Graham et al claim:
“there’s no corresponding chain rule of finite calculus, because there’s no
nice form for $\Delta f(g(x))$.” To the contrary, since $\Delta x=(x+1)-x=1$,
an equivalent way to write $\Delta f(x)$ is

$\Delta f(x)=f(x+\Delta x)-f(x).$

One can then easily check that

$\Delta f(g(x))=f(g(x+\Delta x))-f(g(x))=f(g(x)+\Delta g(x))-f(g(x)),$

which can be expressed as a chain rule

$\Delta(f(x)\llbracket x|g(x)\rrbracket)=(\Delta f(x))\llbracket
x|g(x)\rrbracket,$

where $\llbracket x|g(x)\rrbracket$ operates by replacing $x$ by $g(x)$ and
$\Delta x$ by $\Delta g(x)$.

Of course, to make this rigorous, it would be necessary to work in a formal
language so as to carefully track which “$1$”s are “$\Delta x$”s. For example,
if $f(x)=1/(1+x^{2})$, we want $f(x)\llbracket x|g(x)\rrbracket$ to be
$1/(1+g(x)^{2})$, not $\Delta g(x)/(\Delta g(x)+g(x)^{2})$, even though
$\Delta x=1$. We will not go through the necessary formalism in this paper,
but it would be very similar to the formalism required for the $d$ operator,
which we devote the whole rest of the paper to.

##  5\. Formalizing terms

In this section, we will formalize the terms (or expressions) of differential
calculus. We attempt to make this formalization self-contained. The machinery
we develop here is very similar to the machinery used to define terms in
first-order logic, except that we assume more structure on the set of
variables than is assumed in first-order logic.

Note that one could strongly argue that elementary calculus already implicitly
operates on terms, abusing language to call terms “functions”. For example,
$x\mapsto x^{2}$ and $y\mapsto y^{2}$ are two names for the exact same
function. Yet, nevertheless, in elementary calculus, the expressions $x^{2}$
and $y^{2}$ are _not_ interchangeable [5]. Evidently, such discrepancies point
to the fact that elementary calculus really is done using formal terms,
implicitly. In the following, we make it explicit.

######  Definition 5.1.

(Variables) We fix a set of _variables_ defined inductively as follows.

  1. (1)

For the base step, we fix a countably infinite set $\\{x_{0},x_{1},\ldots\\}$
of distinct elements called _precalculus variables_ , and we declare them to
be variables.

  2. (2)

Inductively, for every variable $v$, we fix a new variable $dv$, which we call
a _differential variable_ ; we do this in such a way as to satisfy the
following requirement (we write $d^{n}v$ for $ddd\cdots dv$ where $d$ occurs
$n$ times):

     * •

(Unique Readability) For all $n,m\in\mathbb{N}$, for all variables $v$ and
$w$, if $d^{n}v$ is the same variable as $d^{m}w$, then $n=m$ and $v=w$.

We write $\mathscr{V}$ for the set of variables.

Examples of variables include $x_{1}$, $x_{50}$, $dx_{0}$, $ddx_{3}$,
$d^{4}x_{50}$ (shorthand for $ddddx_{50}$), and so on. The unique readability
property guarantees that, for example, $dx_{1}$ is not the same variable as
$dx_{2}$ or $ddx_{3}$ or $dddx_{1}$, etc. We allow $n$ or $m$ to be $0$ in the
unique readability requirement, so, for example, $x_{1}$ and $dx_{1}$ are not
the same variable (since $d^{0}x_{1}$ denotes $x_{1}$). Every variable is
either a precalculus variable (in which case it is $x_{n}$ for some
$n\in\mathbb{N}$) or a differential variable (in which case it is $d^{m}x_{n}$
for some $n,m\in\mathbb{N}$ with $m>0$).

######  Definition 5.2.

(Constant symbols and function symbols)

  1. (1)

We fix a distinct set $\\{\overline{r}\\}_{r\in\mathbb{R}}$ of _constant
symbols_ for the real numbers. For any $r\in\mathbb{R}$, $\overline{r}$ is the
_constant symbol for $r$_.

  2. (2)

For every $n\in\mathbb{N}$ with $n>0$, we fix a distinct set
$\\{\overline{f}\\}_{f}$ of _$n$ -ary function symbols_, where $f$ ranges over
the set of all functions from $\mathbb{R}^{n}$ to $\mathbb{R}$. For any such
$f$, $\overline{f}$ is the _$n$ -ary function symbol for $f$_.

We make these choices in such a way that no variable is a constant symbol, no
variable is an $n$-ary function symbol (for any $n$), and no constant symbol
is an $n$-ary function symbol (for any $n$).

For example, the exponential function $\exp$ gives rise to a $1$-ary (or
_unary_) function symbol $\overline{\exp}$. The addition function $+$ gives
rise to a $2$-ary (or _binary_) function symbol $\overline{+}$.

######  Definition 5.3.

(Terms) We define the _terms of differential calculus_ (or simply _terms_)
inductively as follows.

  1. (1)

Every variable $v$ is a term.

  2. (2)

Every constant symbol is a term.

  3. (3)

For all $n\in\mathbb{N}$ ($n>0$), for every $f:\mathbb{R}^{n}\to\mathbb{R}$,
for all terms $U_{1},\ldots,U_{n}$, $\overline{f}(U_{1},\ldots,U_{n})$ is a
term.

Examples of terms include $\overline{5}$, $\overline{\pi}$, $x_{1}$, $dx_{2}$,
$\overline{\sin}(x_{1})$, $\overline{+}(x_{0},x_{1})$, and so on. We often
abuse notation and suppress the overlines and possibly parentheses when
writing terms. For example, we might write $\sin x_{0}$ instead of
$\overline{\sin}(x_{0})$, $\cos\pi$ instead of
$\overline{\cos}(\overline{\pi})$, and so on. For certain well-known
functions, we sometimes abuse notation further, for example, writing:

  * •

$x_{0}+x_{1}$ instead of $\overline{+}(x_{0},x_{1})$;

  * •

$2x_{0}$ instead of $\overline{\cdot}(\overline{2},x_{0})$;

  * •

$x_{0}\,dx_{1}$ instead of $\overline{\cdot}(x_{0},dx_{1})$;

  * •

$x^{2}_{0}$ instead of $\overline{x\mapsto x^{2}}(x_{0})$;

  * •

$e^{x_{1}}$ instead of $\overline{exp}(x_{1})$;

  * •

$x_{0}\,dx_{1}+x_{1}\,dx_{0}$ instead of
$\overline{+}(\overline{\cdot}(x_{0},dx_{1}),\overline{\cdot}(x_{1},dx_{0}))$;

  * •

and so on.

This should cause no confusion in practice.

######  Definition 5.4.

(Term interpretation)

  * •

By an _assignment_ , we mean a function $s:\mathscr{V}\to\mathbb{R}$ (recall
that $\mathscr{V}$ is the set of variables).

  * •

Let $s$ be an assignment. For every term $T$, we define the _interpretation_
$T^{s}\in\mathbb{R}$ of $T$ (according to $s$) by induction on term complexity
as follows.

    1. (1)

If $T$ is a constant symbol $\overline{r}$, then $T^{s}=r$.

    2. (2)

If $T$ is a variable $v$, then $T^{s}=s(v)$.

    3. (3)

If $T$ is $\overline{f}(U_{1},\ldots,U_{n})$ for some
$f:\mathbb{R}^{n}\to\mathbb{R}$ and terms $U_{1},\ldots,U_{n}$, then
$T^{s}=f(U^{s}_{1},\ldots,U^{s}_{n})$.

For example, if $s(x_{0})=5$, then $\overline{\exp}(x_{0})^{s}=e^{5}$. If
$s(x_{0})=9$ and $s(dx_{0})=0.1$, then $(x\,dx)^{s}=9\cdot 0.1=0.9$.

######  Definition 5.5.

(Free variables) We define the _free variables_ $\mathrm{FV}(T)$ of a term $T$
as follows.

  1. (1)

If $T$ is a constant symbol, then $\mathrm{FV}(T)=\emptyset$ (the empty set).

  2. (2)

If $T$ is a variable $v$, then $\mathrm{FV}(T)=\\{v\\}$.

  3. (3)

If $T$ is $\overline{f}(U_{1},\ldots,U_{n})$ for some
$f:\mathbb{R}^{n}\to\mathbb{R}$ and terms $U_{1},\ldots,U_{n}$, then

$\mathrm{FV}(T)=\mathrm{FV}(U_{1})\cup\cdots\cup\mathrm{FV}(U_{n}).$

For example, $\mathrm{FV}(\overline{5})=\emptyset$,
$\mathrm{FV}(x_{6})=\\{x_{6}\\}$, $\mathrm{FV}(dx_{2})=\\{dx_{2}\\}$ (note
that $x_{2}$ is _not_ a free variable of $dx_{2}$),
$\mathrm{FV}(e^{x_{0}+x_{1}})=\\{x_{0},x_{1}\\}$,
$\mathrm{FV}(x_{1}\,dx_{2})=\\{x_{1},dx_{2}\\}$.

######  Lemma 5.6.

Suppose $T$ is a term, $v$ is a variable, and $s$ is an assignment. If
$v\not\in\mathrm{FV}(T)$, then $T^{s}$ does not depend on $s(v)$.

###### Proof.

By induction. ∎

######  Definition 5.7.

(Semantic equivalence) If $T$ and $U$ are terms, we declare $T\equiv U$ (and
say that $T$ and $U$ are _semantically equivalent_) if for every assignment
$s$, $T^{s}=U^{s}$.

For example, $\sin(x_{0}+2\pi)\equiv\sin x_{0}$, by which we mean
$\overline{\sin}(\overline{+}(x_{0},\overline{2\pi})))\equiv\overline{\sin}(x_{0})$.

###  5.1. Formal derivatives

######  Definition 5.8.

(Ordered free variables) If $T$ is a term, we define the _ordered free
variables_ $\mathrm{OFV}(T)$ to be the finite sequence whose elements are the
free variables $\mathrm{FV}(T)$ of $T$ (each appearing exactly one time in the
sequence), ordered such that:

  * •

Whenever $0<n<m$ then $d^{n}x_{i}$ precedes $d^{m}x_{j}$.

  * •

Whenever $0<i<j$ then $d^{n}x_{i}$ precedes $d^{n}x_{j}$.

For example,

$\mathrm{OFV}(e^{x_{1}+x_{3}+x_{2}+x_{2}+x_{99}}\,dx_{1}\,d^{3}x_{1}\,dx_{2}\,d^{50}x_{0})=(x_{1},x_{2},x_{3},x_{99},dx_{1},dx_{2},d^{3}x_{1},d^{50}x_{0}).$

######  Definition 5.9.

If $s$ is an assignment, $w$ is a variable, and $r\in\mathbb{R}$, we write
$s(w|r)$ for the assignment defined by

$s(w|r)(v)=\begin{cases}r&\mbox{if $v$ is $w$}\\\
s(v)&\mbox{otherwise}.\end{cases}$

In other words, $s(w|r)$ is the assignment which is identical to $s$ except
that it overrides $s$’s output on $w$, mapping $w$ to $r$ instead.

######  Lemma 5.10.

For any assignment $s$ and variable $v$, $s(v|s(v))=s$.

###### Proof.

Trivial. ∎

######  Definition 5.11.

(Everywhere-differentiability) Let $T$ be a term, $w$ a variable. We say that
_$T$ is everywhere-differentiable with respect to $w$_ if for every assignment
$s$, the limit

$\lim_{h\to 0}\frac{T^{s(w|s(w)+h)}-T^{s}}{h}$

converges to a finite real number.

######  Lemma 5.12.

Let $T$ be a term with $\mathrm{OFV}(T)=(v_{1},\ldots,v_{n})\not=\emptyset$,
and let $w$ be a variable. Assume $T$ is everywhere-differentiable with
respect to $w$. For all $r_{1},\ldots,r_{n}$, let

$f(r_{1},\ldots,r_{n})=\lim_{h\to 0}\frac{T^{s(w|s(w)+h)}-T^{s}}{h}$

where $s$ is some assignment such that each $s(v_{i})=r_{i}$. Then
$f:\mathbb{R}^{n}\to\mathbb{R}$ is well-defined.

###### Proof.

In other words, for any $r_{1},\ldots,r_{n}\in\mathbb{R}$,
$f(r_{1},\ldots,r_{n})$ does not depend on the choice of $s$, as long as each
$s(v_{i})=r_{i}$. This follows from Lemma 5.6 since $T$ has no free variables
other than $v_{1},\ldots,v_{n}$. ∎

######  Definition 5.13.

If $T$ is a term with $\mathrm{OFV}(T)=(v_{1},\ldots,v_{n})$, $w$ is a
variable, and $T$ is everywhere-differentiable with respect to $w$, then we
define the _derivative of $T$ with respect to $w$_, a term, written
$\frac{\partial T}{\partial w}$, as

$\frac{\partial T}{\partial w}=\overline{f}(v_{1},\ldots,v_{n})$

where $f$ is as in Lemma 5.12. We define $\frac{\partial T}{\partial w}$ to be
the term $\overline{0}$ if $\mathrm{FV}(T)=\emptyset$.

######  Example 5.14.

(Some example term derivatives)

  1. (1)

$\partial x_{0}/\partial x_{0}\equiv 1$.

  2. (2)

$\partial x_{0}/\partial x_{1}\equiv 0$.

  3. (3)

$\partial x_{0}/\partial dx_{0}\equiv 0$.

  4. (4)

$\partial(e^{x_{1}x_{2}}\,dx_{1})/\partial x_{1}\equiv
x_{2}e^{x_{1}x_{2}}\,dx_{1}$.

###### Proof.

(1) The function $f$ of Lemma 5.12 is

$f(r)=\lim_{h\to 0}\frac{x_{0}^{s(x_{0}|s(x_{0})+h)}-x_{0}^{s}}{h}$

(for any assignment $s$ with $s(x_{0})=r$). By Definitions 5.4 and 5.9 this
simplifies to $f(r)=\lim_{h\to 0}\frac{s(x_{0})+h-s(x_{0})}{h}=1$. The claim
follows.

(2) The function $f$ of Lemma 5.12 is

$f(r)=\lim_{h\to 0}\frac{x_{0}^{s(x_{1}|s(x_{1})+h)}-x_{0}^{s}}{h}$

(where $s(x_{0})=r$). This simplifies to $f(r)=\lim_{h\to
0}\frac{s(x_{0})-s(x_{0})}{h}=0$. The claim follows.

(3) Similar to (2).

(4) By Definition 5.8,
$\mathrm{OFV}(e^{x_{1}x_{2}}\,dx_{1})=(x_{1},x_{2},dx_{1})$. So, letting
$v_{1}=x_{1}$, $v_{2}=x_{2}$, $v_{3}=dx_{1}$, the function $f$ of Definition
5.12 is

$f(r_{1},r_{2},r_{3})=\lim_{h\to
0}\frac{(e^{x_{1}x_{2}}\,dx_{1})^{s(v_{1}|s(v_{1})+h)}-(e^{x_{1}x_{2}}\,dx_{1})^{s}}{h}$

(where each $s(v_{i})=r_{i}$). By Definitions 5.4 and 5.9 this simplifies to

$f(r_{1},r_{2},r_{3})=\lim_{h\to
0}\frac{e^{(r_{1}+h)r_{2}}r_{3}-e^{r_{1}r_{2}}r_{3}}{h},$

which is $r_{2}e^{r_{1}r_{2}}r_{3}$ by calculus. The claim follows. ∎

Another way to prove Example 5.14 would be to use the following lemma.

######  Lemma 5.15.

For each term $T$, variable $w$, and assignment $t$, if $T$ is everywhere-
differentiable with respect to $w$, then

$\left(\frac{\partial T}{\partial w}\right)^{t}=\lim_{h\to
0}\frac{T^{t(w|t(w)+h)}-T^{t}}{h}.$

###### Proof.

If $\mathrm{FV}(T)=\emptyset$, the lemma is trivial. Assume not. Let
$(v_{1},\ldots,v_{n})=\mathrm{OFV}(T)$. By definition, $\frac{\partial
T}{\partial w}=\overline{f}(v_{1},\ldots,v_{n})$, where
$f:\mathbb{R}^{n}\to\mathbb{R}$ is such that for all
$r_{1},\ldots,r_{n}\in\mathbb{R}$, for any assignment $s$ with each
$s(v_{i})=r_{i}$,

$f(r_{1},\ldots,r_{n})=\lim_{h\to 0}\frac{T^{s(w|s(w)+h)}-T^{s}}{h}.$

In particular, let each $r_{i}=t(v_{i})$. Then:

$\displaystyle\left(\mbox{$\frac{\partial T}{\partial w}$}\right)^{t}$
$\displaystyle=\overline{f}(v_{1},\ldots,v_{n})^{t}$ (Definition 5.13)
$\displaystyle=f(t(v_{1}),\ldots,t(v_{n}))$ (Definition 5.4)
$\displaystyle=f(r_{1},\ldots,r_{n})$ (Choice of $r_{1},\ldots,r_{n}$)
$\displaystyle=\lim_{h\to 0}\frac{T^{t(w|t(w)+h)}-T^{t}}{h},$ (Since each
$t(v_{i})=r_{i}$)

as desired. ∎

######  Definition 5.16.

(Term total differentials) Suppose $T$ is a term. We say $T$ is _everywhere
totally differentiable_ if $T$ is everywhere-differentiable with respect to
every variable. If so, we define the _total differential_ $\mathbf{d}T$, a
term, as follows. If $\mathrm{FV}(T)=\emptyset$ then we define
$\mathbf{d}T=\overline{0}$. Otherwise, let
$\mathrm{OFV}(T)=(v_{1},\ldots,v_{n})$ and define

$\mathbf{d}T=\frac{\partial T}{\partial v_{1}}dv_{1}+\cdots+\frac{\partial
T}{\partial v_{n}}dv_{n}.$

Furthermore, we inductively define $\mathbf{d}^{1}T$ to be $\mathbf{d}T$ and,
whenever $\mathbf{d}^{n}T$ is defined and is everywhere totally
differentiable, we define $\mathbf{d}^{n+1}T=\mathbf{d}\mathbf{d}^{n}T$.

For example,

$\displaystyle\mathbf{d}(x_{1}\,dx_{2})$
$\displaystyle=\frac{\partial(x_{1}\,dx_{2})}{\partial
x_{1}}dx_{1}+\frac{\partial(x_{1}\,dx_{2})}{\partial dx_{2}}ddx_{2}$
$\displaystyle\equiv dx_{1}\,dx_{2}+x_{1}\,ddx_{2}.$

######  Lemma 5.17.

If term $T$ is everywhere totally differentiable and if $v_{1},\ldots,v_{n}$
are distinct variables such that
$\mathrm{FV}(T)\subseteq\\{v_{1},\ldots,v_{n}\\}$, then

$\mathbf{d}T\equiv\frac{\partial T}{\partial
v_{1}}dv_{1}+\cdots+\frac{\partial T}{\partial v_{n}}dv_{n}.$

###### Proof.

Follows from the commutativity of addition and the fact that clearly
$\frac{\partial T}{\partial v_{i}}\equiv\overline{0}$ if
$v_{i}\not\in\mathrm{FV}(T)$. ∎

In order to prove an abstract chain rule in Section 6, we will need a form of
the classical multivariable chain rule, expressed for formal terms. For this
purpose, we first introduce shorthand for finite summation notation222It is
also possible to incorporate summation notation formally into Definition 5.3,
but the details are complicated. See [1]..

######  Definition 5.18.

If $m>0$ is an integer and $T_{1},\ldots,T_{m}$ are terms, we write
$\sum_{i=1}^{m}T_{i}$ (or just $\sum_{i}T_{i}$ if no confusion results) as
shorthand for $T_{1}+\cdots+T_{m}$.

######  Lemma 5.19.

(Classic Multivariable Chain Rule for Terms) Suppose
$f:\mathbb{R}^{n}\to\mathbb{R}$. Suppose $\vec{T}=(T_{1},\ldots,T_{n})$ are
terms with each $\mathrm{FV}(T_{i})\subseteq\\{v_{1},\ldots,v_{m}\\}$ (where
$v_{1},\ldots,v_{m}$ are distinct). Assume that $\overline{f}(\vec{T})$ and
$T_{1},\ldots,T_{n}$ are everywhere totally differentiable. Then for all
$1\leq i\leq m$,

$\frac{\partial(\overline{f}(\vec{T}))}{\partial
v_{i}}\equiv\sum_{j=1}^{n}\overline{f_{j}}(\vec{T})\frac{\partial
T_{j}}{\partial v_{i}},$

where $f_{j}=D_{j}f$ (the partial derivative of $f$ (in the usual sense) with
respect to its $j$th argument).

###### Proof.

Let $s$ be an assignment and fix $1\leq i\leq m$. We must show (Definition
5.7) that

$\left(\frac{\partial(\overline{f}(\vec{T}))}{\partial
v_{i}}\right)^{s}=\left(\sum_{j=1}^{n}\overline{f_{j}}(\vec{T})\frac{\partial
T_{j}}{\partial v_{i}}\right)^{s}.$

Define functions $F,G_{j}:\mathbb{R}\to\mathbb{R}$ ($1\leq j\leq n$) by

$\displaystyle F(z)$ $\displaystyle=\overline{f}(\vec{T})^{s(v_{i}|z)},$
$\displaystyle G_{j}(z)$ $\displaystyle=T^{s(v_{i}|z)}_{j}.$

For all $1\leq j\leq n$ and $z\in\mathbb{R}$,

$\displaystyle F(z)$ $\displaystyle=\overline{f}(\vec{T})^{s(v_{i}|z)}$
(Definition of $F_{i}$)
$\displaystyle=f(T^{s(v_{i}|z)}_{1},\ldots,T^{s(v_{i}|z)}_{n})$ (Definition
5.4) $\displaystyle=f(G_{1}(z),\ldots,G_{n}(z)),$ (Definition of $G_{j}$)
$\displaystyle\mbox{so ($*$) }F^{\prime}(z)$
$\displaystyle=\mbox{$\sum_{j}$}f_{j}(G_{1}(z),\ldots,G_{n}(z))G^{\prime}_{j}(z)$
(Classic multivar. chain rule)

(the hypotheses of the classic multivariable chain rule are implied by the
everywhere-total-differentiability of $\overline{f}(\vec{T})$ and each
$T_{i}$, by Lemma 5.15). So armed, we compute:

$\displaystyle\left(\mbox{$\frac{\partial(\overline{f}(\vec{T}))}{\partial
v_{i}}$}\right)^{s}$ $\displaystyle=\lim_{h\to
0}\frac{\overline{f}(\vec{T})^{s(v_{i}|s(v_{i})+h)}-\overline{f}(\vec{T})^{s}}{h}$
(Lemma 5.15) $\displaystyle=\lim_{h\to 0}\frac{F(s(v_{i})+h)-F(s(v_{i}))}{h}$
(Def. of $F$) $\displaystyle=F^{\prime}(s(v_{i}))$ (Def. of $F^{\prime}$)
$\displaystyle=\mbox{$\sum_{j}$}f_{j}(G_{1}(s(v_{i})),\ldots,G_{n}(s(v_{i})))G^{\prime}_{j}(s(v_{i}))$
(By ($*$))
$\displaystyle=\mbox{$\sum_{j}$}f_{j}(T^{s(v_{i}|s(v_{i}))}_{1},\ldots,T^{s(v_{i}|s(v_{i}))}_{n})G^{\prime}_{j}(s(v_{i}))$
(Def. of $G_{j}$)
$\displaystyle=\mbox{$\sum_{j}$}f_{j}(T^{s}_{1},\ldots,T^{s}_{n})G^{\prime}_{j}(s(v_{i}))$
(Lemma 5.10)
$\displaystyle=\mbox{$\sum_{j}$}f_{j}(T^{s}_{1},\ldots,T^{s}_{n})\lim_{h\to
0}\frac{G_{j}(s(v_{i})+h)-G_{j}(s(v_{i}))}{h}$ (Def. of $G^{\prime}_{j}$)
$\displaystyle=\mbox{$\sum_{j}$}f_{j}(T^{s}_{1},\ldots,T^{s}_{n})\lim_{h\to
0}\frac{T^{s(v_{i}|s(v_{i})+h)}-T^{s(v_{i}|s(v_{i}))}}{h}$ (Def. of $G_{j}$)
$\displaystyle=\mbox{$\sum_{j}$}f_{j}(T^{s}_{1},\ldots,T^{s}_{n})\lim_{h\to
0}\frac{T^{s(v_{i}|s(v_{i})+h)}-T^{s}}{h}$ (Lemma 5.10)
$\displaystyle=\mbox{$\sum_{j}$}f_{j}(T^{s}_{1},\ldots,T^{s}_{n})\left(\mbox{$\frac{\partial
T_{j}}{\partial v_{i}}$}\right)^{s}$ (Lemma 5.15)
$\displaystyle=\left(\mbox{$\sum_{j=1}^{n}$}\overline{f_{j}}(\vec{T})\mbox{$\frac{\partial
T_{j}}{\partial v_{i}}$}\right)^{s},$ (Def. 5.4)

as desired. ∎

Note that in Lemma 5.19 the assumption that $\overline{f}(\vec{T})$ is
everywhere totally differentiable does not automatically imply that
$T_{1},\ldots,T_{n}$ are everywhere totally differentiable. For example, $f$
could be the function $f(x,y)=x$ in which case $f(T_{1},T_{2})$ would be
everywhere totally differentiable iff $T_{1}$ is everywhere totally
differentiable, regardless of the behavior of $T_{2}$.

##  6\. An Abstract Chain Rule

Recall that $\mathscr{V}$ denotes the set of all variables. Let $\mathscr{T}$
denote the set of all terms.

######  Definition 6.1.

For any $\phi_{0}:\mathscr{V}\to\mathscr{T}$, the _extension of $\phi_{0}$ to
all terms_ is the function $\phi:\mathscr{T}\to\mathscr{T}$ defined by
induction as follows:

  1. (1)

If $T$ is a constant symbol then $\phi(T)=T$.

  2. (2)

If $T$ is a variable then $\phi(T)=\phi_{0}(T)$.

  3. (3)

If $T$ is $\overline{f}(S_{1},\ldots,S_{n})$ then
$\phi(T)=\overline{f}(\phi(S_{1}),\ldots,\phi(S_{n}))$.

######  Lemma 6.2.

Let $\phi_{0}:\mathscr{V}\to\mathscr{T}$ and let $\phi$ be the extension of
$\phi_{0}$ to all terms. Then:

  1. (1)

(The Substitution Lemma) For any assignment $s$, if $\phi(s)$ is the
assignment defined by $\phi(s)(v)=\phi(v)^{s}$, then for every term $T$,
$\phi(T)^{s}=T^{\phi(s)}$.

  2. (2)

For all terms $T$ and $U$, if $T\equiv U$ then $\phi(T)\equiv\phi(U)$.

###### Proof.

(1) By induction on $T$. If $T$ is a constant symbol or variable, the claim is
trivial. Otherwise, $T$ is $\overline{f}(U_{1},\ldots,U_{n})$. Then

$\displaystyle\phi(T)^{s}$
$\displaystyle=\overline{f}(\phi(U_{1}),\ldots,\phi(U_{n}))^{s}$ (Definition
6.1) $\displaystyle=f(\phi(U_{1})^{s},\ldots,\phi(U_{n})^{s})$ (Definition
5.4) $\displaystyle=f(U_{1}^{\phi(s)},\ldots,U_{n}^{\phi(s)})$ (Induction)
$\displaystyle=T^{\phi(s)}.$ (Definition 5.4)

(2) Assume $T\equiv U$. For any assignment $s$, if $\phi(s)$ is as in (1),
then $T^{\phi(s)}=U^{\phi(s)}$ by Definition 5.7. Thus
$\phi(T)^{s}=\phi(U)^{s}$ by (1). By arbitrariness of $s$,
$\phi(T)\equiv\phi(U)$. ∎

######  Definition 6.3.

Say $\phi_{0}:\mathscr{V}\to\mathscr{T}$ _respects $d$_ if for each variable
$v$, $\phi_{0}(dv)\equiv\mathbf{d}\phi_{0}(v)$.

######  Definition 6.4.

(Strong differentiability)

  1. (1)

We define the _subterms_ of a term $T$ by induction as follows. If $T$ is a
variable or constant symbol, then $T$ is its own lone subterm. If $T$ is
$\overline{f}(U_{1},\ldots,U_{n})$, then the subterms of $T$ are $T$ itself
along with the subterms of each $U_{i}$.

  2. (2)

A term $T$ is _strongly differentiable_ if every subterm of $T$ is everywhere
totally differentiable.

Thus, a term is strongly differentiable if it is built up from pieces which
are everywhere totally differentiable. An example of a term which is
everywhere totally differentiable but not strongly differentiable is
$|x_{0}|^{2}$, which is everywhere totally differentiable despite having a
subterm $|x_{0}|$ which is not. Note that the ordinary chain rule for
$f(g(x))^{\prime}$ fails when $f(x)=x^{2}$ and $g(x)=|x|$ (these functions
fail the chain rule’s hypotheses): $(|x|^{2})^{\prime}=2x$, but $|x|^{\prime}$
is undefined at $x=0$. We avoid such traps in the following theorem by
requiring strong differentiability.

######  Theorem 6.5.

(General Abstract Chain Rule) Let $\phi_{0}:\mathscr{V}\to\mathscr{T}$ and
assume that $\phi_{0}(v)$ is strongly differentiable for every variable $v$.
Let $\phi$ be the extension of $\phi_{0}$ to all terms. If $T$ is strongly
differentiable and $\phi_{0}$ respects $d$, then
$\mathbf{d}\phi(T)\equiv\phi(\mathbf{d}T)$.

###### Proof.

By induction on $T$. If $T$ is a constant symbol, the theorem is trivial. If
$T$ is a variable, the theorem reduces to the statement that $\phi_{0}$
respects $d$, which is one of the hypotheses. It remains to consider the case
when $T$ is $\overline{f}(\vec{T})$ where $f:\mathbb{R}^{m}\to\mathbb{R}$ and
$\vec{T}=T_{1},\ldots,T_{m}$ are simpler terms. Then $T_{1},\ldots,T_{m}$ are
subterms of $T$, so, since $T$ is strongly differentiable, it follows that
$T_{1},\ldots,T_{m}$ are strongly differentiable. By induction, each
$\mathbf{d}\phi(T_{i})\equiv\phi(\mathbf{d}T_{i})$. Let
$\\{v_{1},\ldots,v_{\ell}\\}=\mathrm{FV}(\phi(T_{1}))\cup\cdots\cup\mathrm{FV}(\phi(T_{m}))$.
For the rest of the proof, whenever $S$ is a term and $v$ is a variable, we
will write $S_{v}$ for $\frac{\partial S}{\partial v}$. Let
$\overrightarrow{\phi(T)}$ denote $\phi(T_{1}),\ldots,\phi(T_{m})$. We
calculate:

$\displaystyle\mathbf{d}\phi(\overline{f}(\vec{T}))$
$\displaystyle\equiv\mbox{$\sum_{i=1}^{\ell}$}\phi(\overline{f}(\vec{T}))_{v_{i}}dv_{i}$
(Lemma 5.17)
$\displaystyle=\mbox{$\sum_{i}$}\overline{f}(\overrightarrow{\phi(T)})_{v_{i}}dv_{i}$
(Definition 6.1)
$\displaystyle\equiv\mbox{$\sum_{i}$}\mbox{$\sum_{j=1}^{m}$}\overline{f_{j}}(\overrightarrow{\phi(T)})\phi(T_{j})_{v_{i}}dv_{i}$
(Lemma 5.19)
$\displaystyle\equiv\mbox{$\sum_{j}$}\overline{f_{j}}(\overrightarrow{\phi(T)})\mbox{$\sum_{i}$}\phi(T_{j})_{v_{i}}dv_{i}$
(Basic algebra)
$\displaystyle\equiv\mbox{$\sum_{j}$}\overline{f_{j}}(\overrightarrow{\phi(T)})\mathbf{d}\phi(T_{j})$
(Lemma 5.17)
$\displaystyle\equiv\mbox{$\sum_{j}$}\overline{f_{j}}(\overrightarrow{\phi(T)})\phi(\mathbf{d}T_{j})$
(Induction Hypothesis)
$\displaystyle=\phi\left(\mbox{$\sum_{j}$}\overline{f_{j}}(\vec{T})\mathbf{d}T_{j}\right)$
(Definition 6.1)
$\displaystyle\equiv\phi\left(\mbox{$\sum_{j}$}\overline{f_{j}}(\vec{T})\mbox{$\sum_{i=1}^{\ell}$}(T_{j})_{v_{i}}\,dv_{i}\right)$
(Lemma 5.17)
$\displaystyle\equiv\phi\left(\mbox{$\sum_{i}$}\mbox{$\sum_{j}$}\overline{f_{j}}(\vec{T})(T_{j})_{v_{i}}\,dv_{i}\right)$
(Basic algebra)
$\displaystyle\equiv\phi(\mbox{$\sum_{i}$}\overline{f}(\vec{T})_{v_{i}}dv_{i})$
(Lemma 5.19) $\displaystyle\equiv\phi(\mathbf{d}\overline{f}(\vec{T}))$ (Lemma
5.17)

(in the last few lines, we use Lemma 6.2 part 2). ∎

A weakness of the familiar chain rule is that it does not iterate. The
following corollary shows that the abstract chain rule does iterate.

######  Corollary 6.6.

For all $\phi_{0}$, $\phi$ and $T$ as in Theorem 6.5, for all $k\in\mathbb{N}$
($k>0$), if $\mathbf{d}^{\ell}T$ exists and is strongly differentiable for all
$\ell<k$, then

$\mathbf{d}^{k}\phi(T)\equiv\phi(\mathbf{d}^{k}T).$

###### Proof.

By repeated applications of Theorem 6.5. ∎

In Sections 2 and 3 we used a special case of Theorem 6.5 which we will now
formalize. Recall that a precalculus variable is one that is not of the form
$dv$ for any variable $v$.

######  Definition 6.7.

(Variable substitution respecting differentials) Let $v$ be a precalculus
variable, $U$ a term such that $\mathbf{d}^{k}U$ is strongly differentiable
for all $k$. For every term $T$, we will define the _result of substituting
$U$ for $v$ in $T$ while respecting differentials_, written $T[v|U]$, as
follows. First, we define $\phi_{0}:\mathscr{V}\to\mathscr{T}$ so that:

  1. (1)

$\phi_{0}(v)=U$.

  2. (2)

For every $k>0$, $\phi_{0}(d^{k}v)=\mathbf{d}^{k}U$.

  3. (3)

For all variables $w$ not of either of the above two forms, $\phi_{0}(w)=w$.

We define $T[v|U]$ to be $\phi(T)$ where $\phi$ is the extension of $\phi_{0}$
to all terms (Definition 6.1).

######  Corollary 6.8.

(Abstract Chain Rule) Let $U,v$ be as in Definition 6.7. If term $T$ is
strongly differentiable, then

$\mathbf{d}(T[v|U])\equiv(\mathbf{d}T)[v|U].$

###### Proof.

If $\phi_{0}$ is as in Definition 6.7 then evidently $\phi_{0}$ satisfies the
hypotheses of Theorem 6.5. The corollary then immediately follows from Theorem
6.5. ∎

######  Corollary 6.9.

(Iterated Abstract Chain Rule) Let $v,T,U$ be as in Corollary 6.8. For all
$k>0$, if $\mathbf{d}^{\ell}T$ is strongly differentiable for all $\ell<k$,
then

$\mathbf{d}^{k}(T[v|U])\equiv(\mathbf{d}^{k}T)[v|U].$

###### Proof.

By repeated applications of Corollary 6.8. ∎

## Acknowledgments

We gratefully acknowledge Bryan Dawson, Tevian Dray, and the reviewers for
generous comments and feedback.

## References

  * [1] Samuel Alexander.  The first-order syntax of variadic functions.  Notre Dame Journal of Formal Logic, 54(1):47–59, 2013.
  * [2] Samuel Alexander.  A one-sentence elementary proof of the combinatorial Faà di Bruno’s formula.  arXiv preprint 2206.02031, 2022.
  * [3] Alex Craik.  Prehistory of Faà di Bruno’s formula.  The American Mathematical Monthly, 112(2):119–130, 2005.
  * [4] Ronald L Graham, Donald E Knuth, and Oren Patashnik.  Concrete Mathematics: A Foundation for Computer Science.  Addison-Wesley, 2nd edition, 1994.
  * [5] Joel David Hamkins.  The differential operator $\frac{d}{dx}$ binds variables.  In Joel David Hamkins: mathematics and philosophy of the infinite (blog). 2012.
  * [6] Michael Hardy.  Combinatorics of partial derivatives.  The Electronic Journal of Combinatorics, 2006.
  * [7] HN Huang, SAM Marcantognini, and NJ Young.  Chain rules for higher derivatives.  The Mathematical Intelligencer, 28(2):61–69, 2006.
  * [8] Warren P Johnson.  The curious history of Faà di Bruno’s formula.  The American Mathematical Monthly, 109(3):217–234, 2002.
  * [9] Tsoy-Wo Ma.  Higher chain formula proved by combinatorics.  The Electronic Journal of Combinatorics, 2009.
  * [10] Christina Osborne and Amelia Tebbe.  A first step toward higher order chain rules in abelian functor calculus.  In Association for Women in Mathematics Research Symposium, pages 97–119. Springer, 2017.
  * [11] James Stewart.  Calculus.  Brooks Cole, 8th edition, 2015.


* (79) P. Lucey, J. F. Cohn, K. M. Prkachin, P. E. Solomon, and I. Matthews, “Painful data: The UNBC-McMaster shoulder pain expression archive database,” in _2011 IEEE International Conference on Automatic Face and Gesture Recognition (FG)_ , 2011, pp. 57–64.
  * (80) F. Sattler, S. Wiedemann, K.-R. Müller, and W. Samek, “Robust and communication-efficient federated learning from non-IID data,” _IEEE transactions on neural networks and learning systems_ , vol. 31, no. 9, pp. 3400–3413, 2020.
  * (81) T. S. Brisimi, R. Chen, T. Mela, A. Olshevsky, I. C. Paschalidis, and W. Shi, “Federated learning of predictive models from federated electronic health records,” _International journal of medical informatics_ , vol. 112, pp. 59–67, 2018.
  * (82) K. Wei, J. Li, M. Ding, C. Ma, H. H. Yang, F. Farokhi, S. Jin, T. Q. S. Quek, and H. V. Poor, “Federated Learning With Differential Privacy: Algorithms and Performance Analysis,” _IEEE Transactions on Information Forensics and Security_ , vol. 15, pp. 3454–3469, 2020.
  * (83) D. Liu, T. Miller, R. Sayeed, and K. D. Mandl, “FADL: Federated-autonomous deep learning for distributed electronic health record,” _arXiv preprint arXiv:1811.11400_ , 2018.
  * (84) P. Vepakomma, O. Gupta, T. Swedish, and R. Raskar, “Split learning for health: Distributed deep learning without sharing raw patient data,” _arXiv preprint arXiv:1812.00564_ , 2018.
  * (85) R. Shao, H. He, H. Liu, and D. Liu, “Stochastic channel-based federated learning for medical data privacy preserving,” _arXiv preprint arXiv:1910.11160_ , 2019.
  * (86) L. Huang, A. L. Shea, H. Qian, A. Masurkar, H. Deng, and D. Liu, “Patient clustering improves efficiency of federated machine learning to predict mortality and hospital stay time using distributed electronic medical records,” _Journal of biomedical informatics_ , vol. 99, p. 103291, 2019\.
  * (87) L. Huang, Y. Yin, Z. Fu, S. Zhang, H. Deng, and D. Liu, “LoAdaBoost: Loss-based adaboost federated machine learning with reduced computational complexity on iid and non-iid intensive care data,” _Plos one_ , vol. 15, no. 4, p. e0230706, 2020.
  * (88) S. Boughorbel, F. Jarray, N. Venugopal, S. Moosa, H. Elhadi, and M. Makhlouf, “Federated uncertainty-aware learning for distributed hospital EHR data,” _arXiv preprint arXiv:1910.12191_ , 2019.
  * (89) P. Papadopoulos, W. Abramson, A. J. Hall, N. Pitropakis, and W. J. Buchanan, “Privacy and trust redefined in federated machine learning,” _Machine Learning and Knowledge Extraction_ , vol. 3, no. 2, pp. 333–356, 2021.
  * (90) S. R. Pfohl, A. M. Dai, and K. Heller, “Federated and differentially private learning for electronic health records,” _arXiv preprint arXiv:1911.05861_ , 2019.
  * (91) D. C. Nguyen, P. N. Pathirana, M. Ding, and A. Seneviratne, “Blockchain and edge computing for decentralized emrs sharing in federated healthcare,” in _GLOBECOM 2020-2020 IEEE Global Communications Conference_. IEEE, 2020, pp. 1–6.
  * (92) Y. Chen, X. Qin, J. Wang, C. Yu, and W. Gao, “Fedhealth: A federated transfer learning framework for wearable healthcare,” _IEEE Intelligent Systems_ , vol. 35, no. 4, pp. 83–93, 2020.
  * (93) D. Liu, D. Dligach, and T. Miller, “Two-stage federated phenotyping and patient representation learning,” _arXiv preprint arXiv:1908.05596_ , 2019\.
  * (94) G. Krishna Gudur and S. K. Perepu, “Federated learning with heterogeneous labels and models for mobile activity monitoring,” _arXiv e-prints_ , pp. arXiv–2012, 2020.
  * (95) D. Liu, T. A. Miller, and K. D. Mandl, “Confederated machine learning on horizontally and vertically separated medical data for large-scale health system intelligence,” _arXiv preprint arXiv:1910.02109_ , 2019.
  * (96) X. Xu, H. Peng, L. Sun, Y. Niu, H. Ma, L. Liu, and L. He, “Federated depression detection from multi-source mobile health data,” _arXiv preprint arXiv:2102.09342_ , 2021.
  * (97) X. Tan, C.-C. H. Chang, and L. Tang, “A tree-based federated learning approach for personalized treatment effect estimation from heterogeneous data sources,” _arXiv preprint arXiv:2103.06261_ , 2021.
  * (98) Z. Yan, J. Wicaksana, Z. Wang, X. Yang, and K.-T. Cheng, “Variation-aware federated learning with multi-source decentralized medical image data,” _IEEE Journal of Biomedical and Health Informatics_ , 2020.
  * (99) W. Li, F. Milletarì, D. Xu, N. Rieke, J. Hancox, W. Zhu, M. Baust, Y. Cheng, S. Ourselin, M. J. Cardoso _et al._ , “Privacy-preserving federated brain tumour segmentation,” in _International Workshop on Machine Learning in Medical Imaging_ , 2019, pp. 133–141.
  * (100) S. Silva, B. A. Gutman, E. Romero, P. M. Thompson, A. Altmann, and M. Lorenzi, “Federated learning in distributed medical databases: Meta-analysis of large-scale subcortical brain data,” in _2019 IEEE 16th international symposium on biomedical imaging (ISBI 2019)_ , 2019, pp. 270–274.
  * (101) M. J. Sheller, B. Edwards, G. A. Reina, J. Martin, S. Pati, A. Kotrotsou, M. Milchenko, W. Xu, D. Marcus, R. R. Colen _et al._ , “Federated learning in medicine: facilitating multi-institutional collaborations without sharing patient data,” _Scientific reports_ , vol. 10, no. 1, pp. 1–12, 2020\.
  * (102) U. C. Srivastava, D. Upadhyay, and V. Sharma, “Intracranial hemorrhage detection using neural network based methods with federated learning,” _arXiv preprint arXiv:2005.08644_ , 2020.
  * (103) M. Malekzadeh, B. Hasircioglu, N. Mital, K. Katarya, M. E. Ozfatura, and D. Gündüz, “Dopamine: Differentially private federated learning on medical data,” _arXiv e-prints_ , pp. arXiv–2101, 2021.
  * (104) F. Ucar and D. Korkmaz, “COVIDiagnosis-Net: deep bayes-squeezenet based diagnosis of the coronavirus disease 2019 (COVID-19) from X-ray images,” _Medical Hypotheses_ , vol. 140, p. 109761, 2020.
  * (105) P. Guo, P. Wang, J. Zhou, S. Jiang, and V. M. Patel, “Multi-institutional collaborations for improving deep learning-based magnetic resonance image reconstruction using federated learning,” _arXiv preprint arXiv:2103.02148_ , 2021.
  * (106) M. Y. Lu, D. Kong, J. Lipkova, R. J. Chen, R. Singh, D. F. Williamsona, T. Y. Chena, and F. Mahmood, “Federated learning for computational pathology on gigapixel whole slide images,” _arXiv preprint arXiv:2009.10190_ , 2020.
  * (107) X. Li, Y. Gu, N. Dvornek, L. H. Staib, P. Ventola, and J. S. Duncan, “Multi-site fMRI analysis using privacy-preserving federated learning and domain adaptation: ABIDE results,” _Medical Image Analysis_ , vol. 65, p. 101765, 2020.
  * (108) H. R. Roth, K. Chang _et al._ , “Federated learning for breast density classification: A real-world implementation,” in _Domain Adaptation and Representation Transfer, and Distributed and Collaborative Learning_. Springer, 2020, pp. 181–191.
  * (109) D. Nguyen, M. Ding, P. N. Pathirana, and A. Seneviratne, “Blockchain and AI-based solutions to combat coronavirus (COVID-19)-like epidemics: A survey,” _IEEE Access_ , vol. 9, pp. 95 730–95 753, 2021.
  * (110) Q.-V. Pham, D. C. Nguyen, T. Huynh-The, W.-J. Hwang, and P. N. Pathirana, “Artificial intelligence (AI) and big data for coronavirus (COVID-19) pandemic: A survey on the state-of-the-arts,” _IEEE Access_ , vol. 8, pp. 130 820–130 839, 2020.
  * (111) M. Loey, F. Smarandache, and N. E. M. Khalifa, “Within the lack of chest COVID-19 X-ray dataset: A novel detection model based on GAN and Deep Transfer Learning,” _Symmetry_ , vol. 12, no. 4, p. 651, 2020.
  * (112) R. Kumar, A. A. Khan, S. Zhang, J. Kumar, T. Yang, N. A. Golalirz, Zakria, I. Ali, S. Shafiq, and W. Wang, “Blockchain-federated-learning and deep learning models for COVID-19 detection using CT imaging,” _arXiv:2007.06537 [cs, eess]_ , 2020, arXiv: 2007.06537.
  * (113) B. Liu, B. Yan, Y. Zhou, Y. Yang, and Y. Zhang, “Experiments of Federated Learning for COVID-19 Chest X-ray images,” _arXiv:2007.05592 [cs, eess]_ , 2020, arXiv: 2007.05592.
  * (114) W. Zhang, T. Zhou, Q. Lu, X. Wang, C. Zhu, H. Sun, Z. Wang, S. K. Lo, and F.-Y. Wang, “Dynamic fusion-based Federated Learning for COVID-19 detection,” _IEEE Internet of Things Journal_ , pp. 1–1, 2021.
  * (115) I. Feki, S. Ammar, Y. Kessentini, and K. Muhammad, “Federated learning for COVID-19 screening from Chest X-ray images,” _Applied Soft Computing_ , vol. 106, p. 107330, 2021.
  * (116) F. Qian and A. Zhang, “The value of federated learning during and post-COVID-19,” _International Journal for Quality in Health Care_ , vol. 33, no. mzab010, 2021.
  * (117) Q. Dou, T. Y. So, M. Jiang, Q. Liu, V. Vardhanabhuti, G. Kaissis, Z. Li, W. Si, H. H. C. Lee, K. Yu, Z. Feng, L. Dong, E. Burian, F. Jungmann, R. Braren, M. Makowski, B. Kainz, D. Rueckert, B. Glocker, S. C. H. Yu, and P. A. Heng, “Federated deep learning for detecting COVID-19 lung abnormalities in CT: a privacy-preserving multinational validation study,” _npj Digital Medicine_ , vol. 4, no. 1, p. 60, 2021.
  * (118) D. Yang, Z. Xu, W. Li, A. Myronenko, H. R. Roth, S. Harmon, S. Xu, B. Turkbey, E. Turkbey, X. Wang, W. Zhu, G. Carrafiello, F. Patella, M. Cariati, H. Obinata, H. Mori, K. Tamura, P. An, B. J. Wood, and D. Xu, “Federated semi-supervised learning for COVID region segmentation in chest CT using multi-national data from China, Italy, Japan,” _Medical Image Analysis_ , vol. 70, p. 101992, 2021.
  * (119) “Federated learning for medical imaging,” 2021. [Online]. Available: https://www.intel.com/content/www/us/en/artificial-intelligence/posts/federated-learning-for-medical-imaging.html
  * (120) “Federated learning brings AI with privacy to hospitals,” 2021. [Online]. Available: https://healthcare-in-europe.com/en/news/federated-learning-brings-ai-with-privacy-to-hospitals.html
  * (121) Q. Yang, L. Fan, R. Tong, and A. Lv, “White paper-IEEE federated machine learning,” _IEEE White Paper_ , 2021.
  * (122) C. A. Choquette-Choo, N. Dullerud, A. Dziedzic, Y. Zhang, S. Jha, N. Papernot, and X. Wang, “CaPC learning: Confidential and private collaborative learning,” _arXiv preprint arXiv:2102.05188_ , 2021.
  * (123) C. de Alwis, A. Kalla, Q.-V. Pham, P. Kumar, K. Dev, W.-J. Hwang, and M. Liyanage, “Survey on 6G frontiers: Trends, applications, requirements, technologies and future research,” _IEEE Open Journal of the Communications Society_ , vol. 2, pp. 836–886, 2021.
  * (124) L. Mucchi, S. Jayousi, S. Caputo, E. Paoletti, P. Zoppi, S. Geli, and P. Dioniso, “How 6G technology can change the future wireless healthcare,” in _2020 2nd 6G wireless summit (6G SUMMIT)_. IEEE, 2020, pp. 1–6.
  * (125) R. Hu, Y. Guo, H. Li, Q. Pei, and Y. Gong, “Personalized federated learning with differential privacy,” _IEEE Internet of Things Journal_ , vol. 7, no. 10, pp. 9530–9539, 2020.
  * (126) H. Fereidooni, S. Marchal, M. Miettinen, A. Mirhoseini, H. Möllering, T. D. Nguyen, P. Rieger, A.-R. Sadeghi, T. Schneider, H. Yalame _et al._ , “SAFELearn: secure aggregation for private federated learning,” in _2021 IEEE Security and Privacy Workshops (SPW)_. IEEE, 2021, pp. 56–62.
  * (127) X. Li, M. Jiang, X. Zhang, M. Kamp, and Q. Dou, “Fedbn: Federated learning on non-IID features via local batch normalization,” _arXiv preprint arXiv:2102.07623_ , 2021.
  * (128) Q. Li, Y. Diao, Q. Chen, and B. He, “Federated learning on non-IID data silos: An experimental study,” _arXiv preprint arXiv:2102.02079_ , 2021.


Further author information: (Send correspondence to Ananya Jana)
Ananya Jana: E-mail<EMAIL_ADDRESS>
# Global and Local Interpretation of black-box Machine Learning models to
determine prognostic factors from early COVID-19 data

Ananya Jana  Dept. of Computer Science, Rutgers University, New Jersey, USA
Carlos D. Minacapelli  Dept. of Medicine, Division of Gastroenterology and
Hepatology, Rutgers Robert Wood Johnson Medical School, New Jersey, USA  Vinod
Rustgi  Dept. of Medicine, Division of Gastroenterology and Hepatology,
Rutgers Robert Wood Johnson Medical School, New Jersey, USA  Dimitris Metaxas
Dept. of Computer Science, Rutgers University, New Jersey, USA

###### Abstract

The COVID-19 corona virus has claimed 4.1 million lives, as of July 24, 2021.
A variety of machine learning models have been applied to related data to
predict important factors such as the severity of the disease, infection rate
and discover important prognostic factors. Often the usefulness of the
findings from the use of these techniques is reduced due to lack of method
interpretability. Some recent progress made on the interpretability of machine
learning models has the potential to unravel more insights while using
conventional machine learning models[1, 2, 3]. In this work, we analyze
COVID-19 blood work data with some of the popular machine learning models;
then we employ state-of-the-art post-hoc local interpretability
techniques(e.g.- SHAP, LIME), and global interpretability techniques(e.g. -
symbolic metamodeling) to the trained black-box models to draw interpretable
conclusions. In the gamut of machine learning algorithms, regressions remain
one of the simplest and most explainable models with clear mathematical
formulation. We explore one of the most recent techniques called symbolic
metamodeling to find the mathematical expression of the machine learning
models for COVID-19. We identify Acute Kidney Injury (AKI), initial Albumin
level (ALB_I), Aspartate aminotransferase (AST_I), Total Bilirubin initial
(TBILI) and D-Dimer initial (DIMER) as major prognostic factors of the disease
severity. Our contributions are - (i) uncover the underlying mathematical
expression for the black-box models on COVID-19 severity prediction task (ii)
we are the first to apply symbolic metamodeling to this task, and (iii)
discover important features and feature interactions. Code repository:
https://github.com/ananyajana/interpretable_covid19.

Figure 1: (a)The symbolic metamodeling method takes as input the query data to
the trained black-box model and the query output variable value and performs
symbolic regression to find the best fitting symbolic model expression g(x)
(b) This method transforms the problem of finding out the unknown function in
the model space to finding out the symbolic expression in the symbolic model
space

###### keywords:

COVID-19, Interpretability, Machine Learning, Symbolic metamodel, Mathematical
expression

##  1 INTRODUCTION

COVID-19 or novel corona virus disease is a highly contagious disease spread
through close contact with infected persons. The corona virus pandemic has
resulted in unprecedented challenges to all aspects of our society such as
healthcare, economy and education. A large number of techniques have been
proposed and applied to different COVID-19 related scientific questions such
as identifying the important factors [4, 5, 6], predicting disease severity[7,
8], predicting patient mortality, predicting infection rate[9], the role of
co-morbidities, predicting ICU requirement[10, 11] etc. The recent Deep
Learning models are known to produce high accuracy results for these tasks.
However, most of these models do not offer deep understanding of how the
methods work and do not provide interpretable results.. Yan et. al [12]
propose an interpretable machine learning model - a decision tree with three
features - 1) lactic dehydrogenase (LDH), 2) lymphocyte and 3) high-
sensitivity C-reactive protein (hs-CRP). While this model is very simple, it
does not generalize well, as found out by a later study by Barish et. al[13].
Tsiknakis et. al[14] propose a deep learning network to predict the severity
of the disease from chest X ray images. They interpret the results by applying
attention maps which are further verified by medical experts. Ramchandani et.
al[9] propose a model to forecast the infection rate changes in COVID-19
cases. This work explains the second order feature interactions but lacks the
interpretability at the individual feature level.

Arik et. al[15] propose a compartmental disease modeling technique to forecast
the progression of COVID-19. In this work, the authors explore the evolution
of the different compartments to explain the prediction. Pal et. al[16]
propose an interpretable deep learning models to classify COVID-19 infected
patients cough from non-COVID-19 patients cough. The authors use a feature
embedding and cough embedding module to achieve this goal. Chen et. al[17]
propose a Random Forest based model to classify severe COVID-19 cases from
non-severe ones. They identify a total of 10 features from their 52 feature
dataset. Gemmer et.al[18] propose an alternative fuzzy classifier based
approach for the task of mortality prediction. Matsuyama et.al[19] propose a
convolutional neural network based model for screening with chest CT image
where wavelet coefficients of the entire image are used without cropping.
Zokaeinikoo et. al[20] propose to classify chest radiography images using
convolutional deep neural network and apply an attention mechanism for
interpretability. Tian et. al[21] proposes to use an LSTM and GRU based
architecture for predicting the number of covi19 positive cases in 3-day,
5-day and 7-day window. Their prediction take into account five risk factors
including population size, preventable hospitalization rate, and violent crime
rate. Wu et al.[22] uses local interpretability models like LIME and SHAP and
also visualize the relationship of the identified important features with the
severity, but they also lack a mathematical formulation of the models.

Most of these models focus on local interpretability methods for numerical
models and Class Activation Maps and Saliency Maps for image based machine
learning models. Our main contribution is that we use local interpretability
methods to identify important features contributing to the machine learning
outcome and we also provide a mathematical formulation for the COVID-19
severity models and thereby helping in global interpretability. Our work is
organized as the follows - We first describe the dataset and its preprocessing
and then we inroduce and outline the interpretability techniques. Next we
provide details on the experiments we performed and the results. Finally we
conclude with a brief discussion of the results and their relevance.

Patient Information |  Features
---|---
Basic Information |  Age, Gender, Race/Ethnicity, BMI, Medical History, Prior medications, co-morbidities,
Presenting Symptoms |  Duration of GI Symptoms prior to admission, Co-infections, GI bleed
Method of Diagnosis |  RT-PCR, Qualitative, Isothermal, Non-PCR, Serological Test, Chest computed tomography
Laboratory/Imaging |  WBC, Hgb, Hct, MCV, Platelet, count initial, Neutrophil Abs initial, Lymphocyte Abs initial, Fibrinogen, Initial, D-Dimer, Glucose Initial, BUN initial, Cr initial, K initial, GFR initial, AKI, Na initial, Bicarb initial, Alb initial-worst-day5, Tbili initial, EGD, Initial EKG, Alk Phos-initial-worst-day5, CT chest, ALT initial-worst-day5, CXR, AST initial-worst-day5, LDH inital, CPK initial, Troponin initial, Ferritin Initial, ProBNP, CRP initial, HBa1c, TSH, T4 Free, Lipase initial, ABG, VBG, Urine Blood, Stool studies, IL-6, Autoimmune markers, T-spot, Colonoscopy, Liver biopsy, Abd ultrasound
Treatment |  Mechanical ventilation, SpO2, Duration of GI symptoms after admission, Admitted to ICU, Length of stay in ICU,Died in hospital, Severe Outcomes, Readmission within 30 days
Table 1: Initial Features in the Dataset

##  2 Dataset Description

###  2.1 Data cleaning and Preprocessing

We have a private dataset of 565 patients who were diagnosed with COVID-19 and
underwent treatment at the specific site between January 1, 2020 and April 30,
2020. The data was protected using HIPAA. The data was de-identified before
being shared and is compliant with the IRB guidelines at the relevant
institutes(Robert Wood Johnson hospital and Rutgers University). The data
consists of different types of information such as demographic information,
COVID-19 screening method, information regarding hospital stay, laboratory
findings of different blood constituents, whether the patients had prior co-
morbidities etc as given in the Table 1. We had a total of 92 initial features
in the dataset. Some of these features were nominal variables with multiple
possible values such as Prior Co-morbidities, Prior Medications, symptoms. We
one-hot? encoded these variables. We also combined some of the features when
they did not result in loss of information, e.g. we merged the features Race
and Ethnicity as Race_Ethnicity. For the patients who have been admitted
multiple times, we discarded all the admission records except the
chronologically last one.

###  2.2 Outcome Severity Metric

We define a new severity metric SEVER to measure the severity of the disease
outcome. The severity metric(SEVER) is assigned the value 1 if any of the
following is true: - (i) the patient died in the hospital, (ii) the patient
had been admitted to ICU, (iii) the patient required mechanical ventilation
support. Otherwise the outcome severity is assigned the value 0.

###  2.3 Feature Selection

We had a total of 214 features in the dataset after the initial data cleaning
and one-hot encoding. All of these features might not be important.

Figure 2: Feature Correlation Matrix

Moreover, building a model with all the features could make the model
unnecessarily complex. Hence in this step we discard some of the non-useful
features and select the most important features in the dataset.

We drop the features ADMDATE (date of admission), DISDATE (date of discharge),
T_stay (type of stay), LENG (length of hospital stay), PAT_ID (patient ID),
V_ID (visit ID), ID (another ID), TREATMENT (treatment offered to the
patients) from feature sets. The reason behind excluding these features is
that they are not intrinsic properties/features of a patient.

There are many missing values in our dataset. We exclude the features where
more than 20% of the patients have missing values of the features while
calculating feature importance.

The one-hot encoded binary categorical features with 90% (or more) of the
values from the same category have been discarded due to low variance.

Mutual information is the amount of information about a random variable that
can be obtained by observing another random variable. It is expressed as
$I(X,Y)=\sum\limits_{x,y}p(x,y)log\frac{p(x,y)}{p(x)p(y)}$
Where $p(x,y)$ is the joint distribution and $p(x)$ and $p(y)$ are the
marginal distributions of the two random variables $X$ and $Y$. We then select
the features(attributes) where the mutual information with the outcome
variable is high. There are still missing values in the dataset. We get rid of
those patient records who still have missing features. We normalize of the
values to be in the range [0, 1].

We then calculate the correlation between the features and also the outcome
variable. We discard the features which have very high correlation with any
other feature or the outcome variable. This is to ensure the model is not
unnecessarily a complex model with a large number of features. From the
correlation matrix Fig. 2 we find that Hgb and Hct are highly correlated and
similarly ALT_I and AST_I are highly correlated with correlation value $>0.8$

Variable Name | Medical Meaning | Description
---|---|---
AKI (any time during the hospitalization) | Acute kidney injury | Laboratory/Imaging
ALB_I or Alb initial (3.5 - 5.5 g/dL) | Albumin | Laboratory/Imaging
ALK_I or Alk Phos inItial (45- 115 IU/L) | Alkaline phosphatase | Laboratory/Imaging
AST_I or AST intial (10-55 IU/L) | Aspartate aminotransferase | Laboratory/Imaging
BUN initial (6-23mg/dL) | Blood urea nitrogen | Laboratory/Imaging
CR or Cr initial (.5-1.2mg/dL) | Creatinine | Laboratory/Imaging
CXR2 | bilateral opacities/infiltrates | Laboratory/Imaging
DIMER or D-Dimer initial (0-500 ng/mL) | - | Laboratory/Imaging
EKG or Initial EKG QTc interval number (esp if received Hydroxychloroquine |  |
or Azythromycin) | Electrocardiogram | Laboratory/Imaging
GFR initial ($>$60mL/min/1.73M*M) | Glomerular filtration rate | Laboratory/Imaging
Glucose initial (70-100mg/dL) | - | Laboratory/Imaging
Hct initial (36.7-44.7%) | Hematocrit | Laboratory/Imaging
MCV initial (78.0-99.0fL) | Mean corpuscular volume | Laboratory/Imaging
NA_I or Na initial (136-145mmol/L) | Sodium | Laboratory/Imaging
TBILI or Tbili initial (0.1 - 1.2 mg/dL) | Total Bilirubin | Laboratory/Imaging
WBC initial (4.9 – 10.0 Thousand/ul) | White blood cells | Laboratory/Imaging
RACE_E | Race and Ethnicity | Personal Information
Age | - | Personal Information
Prior Cmorbidities(HTN) | Hypertension | Co-morbidities
SYMP20 or AMS | Altered mental status | Presenting Symptoms
Table 2: Features Selected and their medical meaning

We also notice that SPO and WBC has negative correlation with the outcome
variable SEVER, this is consistent with the medical knowledge as we know that
the oxygen saturation level drops with COVID-19 severity and White Blood
Corpuscle count increases. We discard the feature SPO as SPO is explicitly
indicative of the disease severity and we intend to explore other prognostic
features which may be rather implicit. At this stage we have 392 patients with
20 attributes or features listed in Table. 2

##  3 Local and Global Interpretability Techniques

Local interpretability methods offer explanation of the outcome for individual
data points whereas global interpretability techniques try to explain the
model globally. In this work we use two local interpretability techniques -
SHAP, LIME and one global interpretability technique symbolic metamodelling.
All of these techniques are model agnostic.

###  3.1 Shapley Additive Analysis(SHAP)

This approach[2] is based on the game theoretic concept of Shapley values that
determines how players contributed marginally to the cost and gain. This
method computes the contribution of each feature in the final prediction for
an instance. Each feature is treated as if the player is in a coalition or a
cooperative game. Shapley values tells us how to fairly distribute the payout
or the predicted outcome among the players or features. These values are
computed for each of these features or players. Shapley value is the average
marginal contribution of a feature across all possible coalitions of features.

###  3.2 Local Interpretable Model-agnostic Explanations(LIME)

This method is based on the concept of local surrogate models. A local
surrogate model tries to explain the prediction of individual instances.[3].
In this method, the black-box model is probed multiple times with samples
which are perturbed versions or variations of the samples from the original
dataset. The predictions for these perturbed versions are noted. A modifed
dataset is built with the perturbed samples and the corresponding predictions.
An interpretable model is trained on this modified dataset and this
interpretable model is weighted by the proximity of the sampled instances to
the instances of interest.

###  3.3 Symbolic Metamodeling

Symbolic metamodel[1] can be thought of as a global surrogate model or model
of a model. Metamodels usually offer a white-box or transparent approximation
of the black-box model. The main challenge associated with black-box models is
the lack of knowledge about the underlying function. This method tries to
approximate the underlying function with the help of Meijer-G functions.
Meijer-G functions are a very general class of function which can express most
of the special functions like exponential, logarithm, cosine, hypergeometric
etc. Meijer-G function is defined as a line integral in the complex plane
$\mathcal{L}$ and is dependent on real-valued parameters. It can be expressed
as:

$G_{p,q}^{m,n}\left(\begin{array}[]{c}a_{1},\ldots,a_{p}\\\
b_{1},\ldots,b_{q}\end{array}\middle|z\right)\\\ =\frac{1}{2\pi
i}\int_{\mathcal{L}}z^{s}\frac{\prod_{j=1}^{m}\Gamma(b_{j}-s)\prod_{j=1}^{n}\Gamma(1-a_{j}+s)}{\prod_{j=m+1}^{q}\Gamma(1-b_{j}+s)\prod_{j=n+1}^{p}\Gamma(a_{j}-s)}\,ds$

where $\Gamma$ denote the gamma function. The metamodel proposed in this
work[1] is parameterized by Meijer-G functions which enables the application
of gradient descent algorithms for parameter optimization. This method is a
post-hoc analysis technique. It queries a trained model for the outcome
variable value with respect to a given input data point as seen in Fig. 1. The
query input data point and the model’s outcome value pair is used by the
symbolic metamodel to train itself and perform gradient descent for parameter
optimization. This technique gives a metamodel with a closed form algebraic
solution.

##  4 Experiments on COVID-19 data

###  4.1 Training different machine learning models

We use three different machine learning models for the severity prediction
task - Random Forest, Light Gradient Boosting Mechanism(LightGBM) and Extreme
Gradient Boosting(XGBoost) Trees. We use a 70:30 train-test split on our data
of 392 patient records. The test set contains 118 patient records(82:non-
severe and 36:severe outcomes).

####  4.1.1 Random Forest

Random Forest is an ensemble learning method used for classification and
regression problems. This method constructs a number of decision trees at
training time and the final prediction is the class which is the mode or
mean/average of all those decision trees. Random Forests try to address the
overfitting problem in single Decision Trees. However, as this method produces
a collection of decision trees, it is complicated to interpret the model
although it performs superior to a single decision tree.

####  4.1.2 Light Gradient Boosting Machine or LightGBM

Gradient boosting method is an ensemble method where weak learners are
converted to be better learners over iterations. LightGBM is a high-
performance, decision tree based gradient boosting framework. It splits the
tree leaf-wise according to the best fit at that stage. Since the split is
leaf-wise the loss reduction is usually better than level-wise loss reduction
as in case of other boosting algorithms. This loss is reduced by using
gradient descent.

####  4.1.3 Extreme Gradient Boosting or XGBoost

XGBoost is one of the most widely used gradient boosting method. This method
is dependent on pre-sorted algorithm and histogram based algorithms for
finding the best split and the split is depth-wise which is unlike the level-
wise split in lightGBM. This method is well regularized i.e. it tries to
penalize more complex models by using L1 and L2 loss.

###  4.2 Results

We run our experiments using the three black-box models and note the important
features as shown in Table. 3.

Model | AUC |  Important Features(descending order)
---|---|---
Random Forest | 0.8285 |  ’AKI’, ’WBC’, ’ALB_I’, ’BUN’, ’AST_I’
LightGBM | 0.8333 |  ’AKI’, ’MCV’, ’EKG’, ’ALB_I’, ’CR’,
XGBoost | 0.8196 |  ’AKI’, ’CXR2’, ’Hct’, ’AST_I’, ’Prior_Co1’
Table 3: Important features identified by individual methods

####  4.2.1 Metric

We used Area Under the ROC Curve(AUC) as the metric for measuring the
performance of the black-box models and their white-box approximations.

###  4.3 Interpreting the Models

####  4.3.1 SHAP Analysis

The SHAP analysis of a patient’s outcome prediction via the three different
trained machine learning models is shown in Fig. 3. As the analysis reveals
AKI and ALB_I reduce the outcome value for this particular patient. SHAP
method also gives a global summary based on the instance-wise analysis. In the
SHAP summary plots for the three different methods i.e. Random Forest,
LightGBM and XGBoost in Fig. 4 we see that all of them unanimously put Acute
Kidney Injury as the most important feature to predict disease severity
outcome. We also note that all the summary plots also reveal ALB_I and AST_I
unanimously in the top five major prognostic features. Higher values of AKI,
higher values of AST_I and lower values of ALB_I push the outcome value
higher. This is consistent with medical knowledge[23, 24, 25].

Figure 3: SHAP analysis for a patient with non-severe outcome by different
black-box models(blue:features which reduce the outcome value, pink:features
which try to push the prediction value higher) [zoom in for better view]
Figure 4: SHAP Summary plots for the three black-box models[zoom in for better
view]

####  4.3.2 LIME Analysis

The trained Random Forest model was probed using LIME method, and the local
explanations were generated for the two patients as shown in the Fig. 5. One
of the two patients had a severe outcome and the other one had a non-severe
outcome. The LIME analysis on the Random Forest model also shows that higher
values of Acute Kidney injury(AKI) tries to increase the prediction value i.e.
make the patient outcome severe whereas higher ALB_I value try to make the
patient outcome non-severe. Similarly high AST_I tries to make the patient
outcome severe.

Figure 5: LIME analysis(red:features that decrease the prediction value, and
green:features that increase the prediction value)

####  4.3.3 Symbolic Metamodeling

We provide the symbolic models or mathematical formulations for the three
individual machine learning models in the Table. 4. Here X0,… X19 are the 20
features in our dataset with the mapping as - AKI[X0], ALB_I[X1], ALK_I[X2],
AST_I[X3], Age[X4], BUN[X5], CR[X6], CXR2[X7], DIMER[X8], EKG[X9],
GFR_BIN[X10], Glucose[X11], Hct[X12], MCV[X13], NA_I[X14], Prior_Co1[X15],
RACE_E[X16], SYMP20[X17], TBILI[X18], WBC[X19]. We also delve deeper into the
individual formulas(refer to the code repository) to gain more insights into
the individual black box models as follows:

Model |  Expression(Symbolic Metamodelling) | AUC(Model) | AUC(Symbolic metamodel)
---|---|---|---
Random Forest |  $\frac{1}{5.63e^{0.2293X^{3}_{0}X^{3}_{1}-0.0488X^{3}_{0}X^{3}_{10}+....-0.0007X_{9}}}$ | 0.8285 | 0.7896
LightGBM |  $\frac{1}{2.34e^{0.0209X^{3}_{0}X^{3}_{1}-0.0006X^{3}_{0}X^{3}_{10}+....-0.0001X_{9}}}$ | 0.8333 | 0.7632
XGBoost |  $\frac{1}{14.01e^{0.4132X^{3}_{0}X^{3}_{1}-0.0297X^{3}_{0}X^{3}_{10}+....-.000098X_{9}}}$ | 0.8196 | 0.7517
Table 4: mathematical expressions obtained using symbolic metamodelling

####  4.3.4 Symbolic metamodel for Random forest

The most important five features i.e. values of the coefficients in descending
order are TBILI[X18](0.00919), AST_I[X]3(0.00865), DIMER[X8](0.00312),
EKG[X9](.00071), ALB_I[X2](0.00041) As we can see, AST_I and ALB_I which are
two of the top five features determined by Random Forest(refer Table. 3 are
part of the top five features from the symbolic metamodel as well. The feature
interactions with highest coefficients are SYMP20[X17]AST_I[X3](3.81),
Hct[X12]Age[X4](3.16). AST_I, which is already found to be an important
feature is found to be a part of the most important feature interaction.

####  4.3.5 Symbolic metamodel for lightGBM

The most important five features i.e. values of the coefficients in descending
order are DIMER[X8](0.000555), EKG[X9](0.000188), ALK_I[X2](0.000148),
AST_I[X3](0.000142), TBILI[X18](0.000086) As we can see, EKG which is one of
the top five features determined by LightGBM(refer Table. 3) is part of the
top five features from the its symbolic metamodel as well.

The feature interactions with highest coefficients are
AST_I[X3]DIMER[X8](0.6419), SYMP20[X17]BUN[X5](0.636) We notice the AST_I is
one of the most important features for this metamodel and it is also part of
the most important feature interaction for this metamodel.

####  4.3.6 Symbolic metamodel for XGBoost model

The most important five features i.e. values of the coefficients in descending
order are DIMER[X8](0.06701), AST_I[X]3(0.02893), TBILI[X18](0.00186),
AKI[X0](0.00143), Glucose[X11](0.00121)

As we can see, AST_I and AKI which are two of the top five features determined
by XGBoost(refer Table. 3 are part of the top five features from the its
symbolic metamodel as well. Interestingly one of the important features
decided by all the black-box models i.e. ALB_I, although not present in the
top five features in this metamodel, it is part of one of the two most
important feature interactions - NA_I[X14]WBC[X19(4.71) AKI[X0]ALB_I[X1](4.47)

We also notice that all the symbolic metamodels evaluate TBILI and DIMER to be
important features. This insight is neither captured by the black-box models
nor their local interpretations although medical knowledge confirms this[26,
27]. As we can see from the results in Table. 4, there is some gap between the
performance of the models and the white-box approximations. This is due to the
fact that our dataset is fairly small and usually the search space of Meijer-G
functions is much larger in comparison, however, even with our limited dataset
we could provide a meaningful functional form.

###  4.4 Discussion

In this work we interpreted the COVID-19 outcome severity at individual
patient’s level as well as globally. Our black-box models identified Acute
Kidney Injury (AKI), initial Albumin level (ALB_I) and Aspartate
aminotransferase (AST_I) as important features. Using the local interpretation
methods we explored how these features contributed to the individual outcomes
and finally we provided a global relationship between the features. This
helped us determine which features are important and by how much they
contribute globally to the outcome, exactly how they are related to each other
and their specific interactions. The global formulation also helped us
discover two important features which were neither captured by the black-box
model nor the local interpretability methods - Total Bilirubin initial (TBILI)
and D-Dimer initial( DIMER). Hence the global interpretation significantly
boosted the trust-worthiness of our black-box models. Using the
interpretability method alongwith black-box models helped us identify five
major prognostic features - Acute Kidney Injury (AKI), initial Albumin level
(ALB_I), Aspartate aminotransferase (AST_I), Total Bilirubin initial (TBILI)
and D-Dimer initial( DIMER).

Recently after more covid19 datasets became available privately and publicly,
some works have been done on interpreting the black-box models, our work is
the first to use both global and local interpretability techniques with
machine learning and draw meaningful conclusions from these complementary
methods. All the interpretability techniques used in this work are model
agnostic and can be easily extended to Deep Learning models and get a
mathematical formulation. Given our encouraging results we can conclude that
incorporating interpretability methods in both machine learning and deep
learning methods applied to larger datasets will help us uncover more
interesting insights into this disease.

###  4.5 Acknowledgments

This research was funded based on partial funding to D. Metaxas from NSF:
IIS-1703883, CNS-1747778, CCF-1733843, IIS-1763523, IIS-1849238 – 825536 and
MURI-Z8424104 -440149.

## References

  * [1] Alaa, A. M. and van der Schaar, M., “Demystifying black-box models with symbolic metamodels,” in [Advances in Neural Information Processing Systems ], 11304–11314 (2019).
  * [2] Lundberg, S. M. and Lee, S.-I., “A unified approach to interpreting model predictions,” in [Advances in neural information processing systems ], 4765–4774 (2017).
  * [3] Ribeiro, M. T., Singh, S., and Guestrin, C., “” why should i trust you?” explaining the predictions of any classifier,” in [Proceedings of the 22nd ACM SIGKDD international conference on knowledge discovery and data mining ], 1135–1144 (2016).
  * [4] Schwab, P., Schütte, A. D., Dietz, B., and Bauer, S., “Clinical predictive models for covid-19: Systematic study,” Journal of medical Internet research 22(10), e21439 (2020).
  * [5] Bao, F. S., He, Y., Liu, J., Chen, Y., Li, Q., Zhang, C. R., Han, L., Zhu, B., Ge, Y., Chen, S., et al., “Triaging moderate covid-19 and other viral pneumonias from routine blood tests,” arXiv preprint arXiv:2005.06546 (2020).
  * [6] Ruan, Q., Yang, K., Wang, W., Jiang, L., and Song, J., “Clinical predictors of mortality due to covid-19 based on an analysis of data of 150 patients from wuhan, china,” Intensive care medicine , 1–3 (2020).
  * [7] Liang, W., Liang, H., Ou, L., Chen, B., Chen, A., Li, C., Li, Y., Guan, W., Sang, L., Lu, J., et al., “Development and validation of a clinical risk score to predict the occurrence of critical illness in hospitalized patients with covid-19,” JAMA Internal Medicine (2020).
  * [8] Singer, A. J., Morley, E. J., Meyers, K., Fernandes, R., Rowe, A. L., Viccellio, P., Thode, H. C., Bracey, A., and Henry, M. C., “Cohort of four thousand four hundred four persons under investigation for covid-19 in a new york hospital and predictors of icu care and ventilation,” Annals of emergency medicine (2020).
  * [9] Ramchandani, A., Fan, C., and Mostafavi, A., “Deepcovidnet: An interpretable deep learning model for predictive surveillance of covid-19 using heterogeneous features and their interactions,” IEEE Access 8, 159915–159930 (2020).
  * [10] Chao, H., Fang, X., Zhang, J., Homayounieh, F., Arru, C. D., Digumarthy, S. R., Babaei, R., Mobin, H. K., Mohseni, I., Saba, L., et al., “Integrative analysis for covid-19 patient outcome prediction,” Medical Image Analysis 67, 101844 (2020).
  * [11] Supreeth, P. S., Gabriel, W., Paulina, P., Morgan, C., Laura, N. B., Kathryn, A. H., Crystal, M., Shibani, M., Gregory, R., Yu-Ping, S., et al., “Development and prospective validation of a transparent deep learning algorithm for predicting need for mechanical ventilation,” medRxiv: the preprint server for health sciences (2020).
  * [12] Yan, L., Zhang, H.-T., Goncalves, J., Xiao, Y., Wang, M., Guo, Y., Sun, C., Tang, X., Jing, L., Zhang, M., et al., “An interpretable mortality prediction model for covid-19 patients,” Nature Machine Intelligence , 1–6 (2020).
  * [13] Barish, M., Bolourani, S., Lau, L. F., Shah, S., and Zanos, T. P., “External validation demonstrates limited clinical utility of the interpretable mortality prediction model for patients with covid-19,” Nature Machine Intelligence , 1–3 (2020).
  * [14] Tsiknakis, N., Trivizakis, E., Vassalou, E. E., Papadakis, G. Z., Spandidos, D. A., Tsatsakis, A., Sánchez-García, J., López-González, R., Papanikolaou, N., Karantanas, A. H., et al., “Interpretable artificial intelligence framework for covid-19 screening on chest x-rays,” Experimental and Therapeutic Medicine 20(2), 727–735 (2020).
  * [15] Arik, S., Li, C.-L., Yoon, J., Sinha, R., Epshteyn, A., Le, L., Menon, V., Singh, S., Zhang, L., Nikoltchev, M., et al., “Interpretable sequence learning for covid-19 forecasting,” Advances in Neural Information Processing Systems 33 (2020).
  * [16] Pal, A. and Sankarasubbu, M., “Pay attention to the cough: Early diagnosis of covid-19 using interpretable symptoms embeddings with cough sound signal processing,” arXiv preprint arXiv:2010.02417 (2020).
  * [17] Chen, Y., Ouyang, L., Bao, F. S., Li, Q., Han, L., Zhu, B., Ge, Y., Robinson, P., Xu, M., Liu, J., et al., “An interpretable machine learning framework for accurate severe vs non-severe covid-19 clinical type classification,” Available at SSRN 3638427 (2020).
  * [18] Gatos, I., Tsantis, S., Spiliopoulos, S., Karnabatidis, D., Theotokas, I., Zoumpoulis, P., Loupas, T., Hazle, J. D., and Kagadis, G. C., “Temporal stability assessment in shear wave elasticity images validated by deep learning neural network for chronic liver disease fibrosis stage assessment,” Medical physics 46(5), 2298–2309 (2019).
  * [19] Matsuyama, E. et al., “A deep learning interpretable model for novel coronavirus disease (covid-19) screening with chest ct images,” Journal of Biomedical Science and Engineering 13(07), 140 (2020).
  * [20] Zokaeinikoo, M., Kazemian, P., Mitra, P., and Kumara, S., “Aidcov: An interpretable artificial intelligence model for detection of covid-19 from chest radiography images,” medRxiv (2020).
  * [21] Tian, T., Jiang, Y., Zhang, Y., Li, Z., Wang, X., and Zhang, H., “Covid-net: A deep learning based and interpretable predication model for the county-wise trajectories of covid-19 in the united states,” medRxiv (2020).
  * [22] Wu, H., Ruan, W., Wang, J., Zheng, D., Li, S., Chen, J., Li, K., Chai, X., and Helal, S., “Is ai model interpretable to combat with covid? an empirical study on severity prediction task,” arXiv e-prints , arXiv–2010 (2020).
  * [23] Batlle, D., Soler, M. J., Sparks, M. A., Hiremath, S., South, A. M., Welling, P. A., Swaminathan, S., et al., “Acute kidney injury in covid-19: emerging evidence of a distinct pathophysiology,” Journal of the American Society of Nephrology 31(7), 1380–1383 (2020).
  * [24] Feketea, G. M. and Vlacha, V., “The diagnostic significance of usual biochemical parameters in coronavirus disease 19 (covid-19): Albumin to globulin ratio and crp to albumin ratio,” Frontiers in Medicine 7 (2020).
  * [25] Cai, Q., Huang, D., Yu, H., Zhu, Z., Xia, Z., Su, Y., Li, Z., Zhou, G., Gou, J., Qu, J., et al., “Covid-19: abnormal liver function tests,” Journal of hepatology 73(3), 566–574 (2020).
  * [26] Liu, Z., Li, J., Long, W., Zeng, W., Gao, R., Zeng, G., Chen, D., Wang, S., Li, Q., Hu, D., et al., “Bilirubin levels as potential indicators of disease severity in coronavirus disease patients: a retrospective cohort study,” Frontiers in Medicine 7 (2020).
  * [27] Rostami, M. and Mansouritorghabeh, H., “D-dimer level in covid-19 infection: a systematic review,” Expert review of hematology 13(11), 1265–1275 (2020).