text:
# Combining False Discovery Rates (FDR)?

Suppose we have a population of patients. For each patient, we measure 1000 features and observe if they suffer from any of 10 diseases. We wish to determine which features, if any, are significant in predicting each disease.

Because of the multiple hypotheses (the 1000 features), we need to correct our $$p$$-values. For example, we might focus on the false discover rate and apply the Benjamini-Hochberg procedure

My question is: do we need to treat this as a single multiple hypothesis problem (1000 features $$\times$$ 10 diseases = 10,000 hypotheses), or can we treat it as 10 individual problems (each with 1000 hypotheses)? I'd certainly prefer the latter, since the $$p$$-value correction will suppress fewer terms.

It's pretty clear that to compute the family-wise error rate, we'd need to combine all the hypotheses. I was hoping that the FDR might behave differently, and allow me to analyze each disease independently.

More broadly, I would appreciate any pointers towards different statistical tool in case this general approach is off-base.

UPDATE:

I'm treating Michael Lew's thoughtful answer below as correct, but I subsequently stumbled across a more (statistically) powerful tool for handling FDR. It seems to be very relevant for my problem. Anyone interested in this problem may find these papers helpful:

You are not testing hypotheses, but fishing for interesting findings. There is nothing wrong with that. Do not do things that are often demanded by those who cannot tell the difference between a preliminary scientific investigation and bad statistics. See my commentary on the ASA's statement on P-values for more detail. (It's in the supplementary material and takes a ridiculous amount of clicking to reach, so here is a direct link to a preprint of it I found online.)

'Correction' of P-value for multiplicity burns the power to detect real effects. Never do it unless you have no alternative, where no alternative comes from an inability to do any sort of follow-up and the absence of any corroboration from other data or theory. Do not dichotomise the results into 'significant' and 'not significant', but show all of the observed effect sizes. (I suspect this response will be voted down, but it is not wrong.)

Treat this as a preliminary investigation. Do not adjust the P-values, but rank the statistical 'interestingness' of the features in order of smallness of the P-values. Then follow up with a study designed to investigate just those features that are statistically and/or scientifically interesting. In this case I would say that even if no follow-up study is possible you should publish the raw P-values so that other investigators can use your data as corroboration for their findings.

• Thank you for that link (which eventually lead me to this and several other interesting papers of yours). Your specific suggestion of sharing the raw P-values is definitely a good one. Your broader point about recognizing exploratory statistics is excellent. This question was motivated by a reviewer's request for corrected p-values; instead of providing that, I think I will try to echo some version of your points above. Thank you! (I am left a little mystified as to when one uses the FDR; the "discovery" part lead me to expect this case) Nov 6, 2017 at 1:55
• PS: I updated my question; if you're not already aware of them, you may find the linked papers interesting. Nov 22, 2017 at 16:55


text:
Isotonic regression

An example of isotonic regression

In numerical analysis, isotonic regression (IR) involves finding a weighted least-squares fit $x\in \Bbb{R}^n$ to a vector $a\in \Bbb{R}^n$ with weights vector $w\in \Bbb{R}^n$ subject to a set of non-contradictory constraints of the kind $x_i \ge x_j$.

Such constraints define partial order or total order and can be represented as a directed graph $G=(N,E)$, where N is the set of variables involved, and E is the set of pairs (i, j) for each constraint $x_i \ge x_j$. Thus, the IR problem corresponds to the following quadratic program (QP):

$\min \sum_{i=1}^n w_i (x_i - a_i)^2$ $\text{subject to }x_i\ge x_j~ \text{ for all } (i,j)\in E.$

In the case when $G=(N,E)$ is a total order, a simple iterative algorithm for solving this QP is called the pool adjacent violators algorithm (PAVA). Best and Chakravarti (1990) have studied the problem as an active set identification problem, and have proposed a primal algorithm in O(n), the same complexity as the PAVA, which can be seen as a dual algorithm.[1]

IR has applications in statistical inference, for example, to fit of an isotonic curve to mean experimental results when an order is expected. A benefit of isotonic regression is that it does not assume any form for the target function, such as linearity assumed by linear regression.

Another application is nonmetric multidimensional scaling,[2] where a low-dimensional embedding for data points is sought such that order of distances between points in the embedding matches order of dissimilarity between points. Isotonic regression is used iteratively to fit ideal distances to preserve relative dissimilarity order.

Isotonic regression is also sometimes referred to as monotonic regression. Correctly speaking, isotonic is used when the direction of the trend is strictly increasing, while monotonic could imply a trend that is either strictly increasing or strictly decreasing.

Isotonic regression under the $L_p$ for $p>0$ is defined as follows:

$\min \sum_{i=1}^n w_i |x_i - a_i|^p$ $\mathrm{subject~to~}x_i\ge x_j~ \text{ for all } (i,j)\in E.$

Simply ordered case

To illustrate the above, let $x_1 \leq x_2 \leq \ldots \leq x_n$, and $f(x_1) \leq f(x_2) \leq \ldots \leq f(x_n)$, and $w_i \geq 0$.

The isotonic estimator, $g^*$, minimizes the weighted least squares-like condition:

$\min_g \sum_{i=1}^n w_i (g(x_i) - f(x_i))^2$

Where $g$ is the unknown function we are estimating, and $f$ is a known function.

Software has been developed in the R statistical package for computing isotone (monotonic) regression. [3]

References

1. ^ Best, M.J.; & Chakravarti N. (1990). "Active set algorithms for isotonic regression; a unifying framework". Mathematical Programming 47: 425–439. doi:10.1007/BF01580873.
2. ^ Kruskal, J. B. (1964). "Nonmetric Multidimensional Scaling: A numerical method". Psychometrika 29 (2): 115–129. doi:10.1007/BF02289694.
3. ^ Leeuw, Jan de; Hornik, Kurt; Mair, Patrick (2009). "Isotone Optimization in R: Pool-Adjacent-Violators Algorithm (PAVA) and Active Set Methods". Journal of Statistical Software 32 (5): 1–24. doi:10.18637/jss.v032.i05. ISSN 1548-7660.

format: plain


text:
# cf.CellMethod.equals¶

CellMethod.equals(other, rtol=None, atol=None, verbose=None, ignore_qualifiers=(), ignore_type=False)[source]

Whether two cell method constructs are the same.

Equality is strict by default. This means that:

• the descriptive qualifiers must be the same (see the ignore_qualifiers parameter).

The axes of the cell method constructs are not considered, because they may only be correctly interpreted by the field constructs that contain the cell method constructs in question. They are, however, taken into account when two fields constructs are tested for equality.

Two real numbers x and y are considered equal if |x-y|<=atol+rtol|y|, where atol (the tolerance on absolute differences) and rtol (the tolerance on relative differences) are positive, typically very small numbers. The data type of the numbers is not taken into consideration. See the atol and rtol parameters.

Any type of object may be tested but, in general, equality is only possible with another cell method construct, or a subclass of one. See the ignore_type parameter.

NetCDF elements, such as netCDF variable and dimension names, do not constitute part of the CF data model and so are not checked.

New in version 1.7.0.

Parameters
other:

The object to compare for equality.

atol: float, optional

The tolerance on absolute differences between real numbers. The default value is set by the cfdm.ATOL function.

rtol: float, optional

The tolerance on relative differences between real numbers. The default value is set by the cfdm.RTOL function.

verbose: int or None, optional

If an integer from 0 to 3, corresponding to increasing verbosity (else -1 as a special case of maximal and extreme verbosity), set for the duration of the method call (only) as the minimum severity level cut-off of displayed log messages, regardless of the global configured cfdm.LOG_LEVEL.

Else, if None (the default value), log messages will be filtered out, or otherwise, according to the value of the cfdm.LOG_LEVEL setting.

Overall, the higher a non-negative integer that is set (up to a maximum of 3) the more description that is printed to convey information about differences that lead to inequality.

ignore_qualifiers: sequence of str, optional

The names of qualifiers to omit from the comparison.

ignore_type: bool, optional

Any type of object may be tested but, in general, equality is only possible with another cell method construct, or a subclass of one. If ignore_type is True then CellMethod(source=other) is tested, rather than the other defined by the other parameter.

Returns
bool

Whether the two cell method constructs are equal.

Examples:

>>> c.equals(c)
True
>>> c.equals(c.copy())
True
>>> c.equals('not a cell method')
False


text:
# Computing Integrals using Meijer G-Functions#

This text aims do describe in some detail the steps (and subtleties) involved in using Meijer G-functions for computing definite and indefinite integrals. We shall ignore proofs completely.

## Overview#

The algorithm to compute $$\int f(x) \mathrm{d}x$$ or $$\int_0^\infty f(x) \mathrm{d}x$$ generally consists of three steps:

1. Rewrite the integrand using Meijer G-functions (one or sometimes two).

2. Apply an integration theorem, to get the answer (usually expressed as another G-function).

3. Expand the result in named special functions.

Step (3) is implemented in the function hyperexpand (q.v.). Steps (1) and (2) are described below. Moreover, G-functions are usually branched. Thus our treatment of branched functions is described first.

Some other integrals (e.g. $$\int_{-\infty}^\infty$$) can also be computed by first recasting them into one of the above forms. There is a lot of choice involved here, and the algorithm is heuristic at best.

## Polar Numbers and Branched Functions#

Both Meijer G-Functions and Hypergeometric functions are typically branched (possible branchpoints being $$0$$, $$\pm 1$$, $$\infty$$). This is not very important when e.g. expanding a single hypergeometric function into named special functions, since sorting out the branches can be left to the human user. However this algorithm manipulates and transforms G-functions, and to do this correctly it needs at least some crude understanding of the branchings involved.

To begin, we consider the set $$\mathcal{S} = \{(r, \theta) : r > 0, \theta \in \mathbb{R}\}$$. We have a map $$p: \mathcal{S}: \rightarrow \mathbb{C}-\{0\}, (r, \theta) \mapsto r e^{i \theta}$$. Decreeing this to be a local biholomorphism gives $$\mathcal{S}$$ both a topology and a complex structure. This Riemann Surface is usually referred to as the Riemann Surface of the logarithm, for the following reason: We can define maps $$\operatorname{Exp}: \mathbb{C} \rightarrow \mathcal{S}, (x + i y) \mapsto (\exp(x), y)$$ and $$\operatorname{Log}: \mathcal{S} \rightarrow \mathbb{C}, (e^x, y) \mapsto x + iy$$. These can both be shown to be holomorphic, and are indeed mutual inverses.

We also sometimes formally attach a point “zero” ($$0$$) to $$\mathcal{S}$$ and denote the resulting object $$\mathcal{S}_0$$. Notably there is no complex structure defined near $$0$$. A fundamental system of neighbourhoods is given by $$\{\operatorname{Exp}(z) : \Re(z) < k\}$$, which at least defines a topology. Elements of $$\mathcal{S}_0$$ shall be called polar numbers. We further define functions $$\operatorname{Arg}: \mathcal{S} \rightarrow \mathbb{R}, (r, \theta) \mapsto \theta$$ and $$|.|: \mathcal{S}_0 \rightarrow \mathbb{R}_{>0}, (r, \theta) \mapsto r$$. These have evident meaning and are both continuous everywhere.

Using these maps many operations can be extended from $$\mathbb{C}$$ to $$\mathcal{S}$$. We define $$\operatorname{Exp}(a) \operatorname{Exp}(b) = \operatorname{Exp}(a + b)$$ for $$a, b \in \mathbb{C}$$, also for $$a \in \mathcal{S}$$ and $$b \in \mathbb{C}$$ we define $$a^b = \operatorname{Exp}(b \operatorname{Log}(a))$$. It can be checked easily that using these definitions, many algebraic properties holding for positive reals (e.g. $$(ab)^c = a^c b^c$$) which hold in $$\mathbb{C}$$ only for some numbers (because of branch cuts) hold indeed for all polar numbers.

As one peculiarity it should be mentioned that addition of polar numbers is not usually defined. However, formal sums of polar numbers can be used to express branching behaviour. For example, consider the functions $$F(z) = \sqrt{1 + z}$$ and $$G(a, b) = \sqrt{a + b}$$, where $$a, b, z$$ are polar numbers. The general rule is that functions of a single polar variable are defined in such a way that they are continuous on circles, and agree with the usual definition for positive reals. Thus if $$S(z)$$ denotes the standard branch of the square root function on $$\mathbb{C}$$, we are forced to define

$\begin{split}F(z) = \begin{cases} S(p(z)) &: |z| < 1 \\ S(p(z)) &: -\pi < \operatorname{Arg}(z) + 4\pi n \le \pi \text{ for some } n \in \mathbb{Z} \\ -S(p(z)) &: \text{else} \end{cases}.\end{split}$

(We are omitting $$|z| = 1$$ here, this does not matter for integration.) Finally we define $$G(a, b) = \sqrt{a}F(b/a)$$.

## Representing Branched Functions on the Argand Plane#

Suppose $$f: \mathcal{S} \to \mathbb{C}$$ is a holomorphic function. We wish to define a function $$F$$ on (part of) the complex numbers $$\mathbb{C}$$ that represents $$f$$ as closely as possible. This process is knows as “introducing branch cuts”. In our situation, there is actually a canonical way of doing this (which is adhered to in all of SymPy), as follows: Introduce the “cut complex plane” $$C = \mathbb{C} \setminus \mathbb{R}_{\le 0}$$. Define a function $$l: C \to \mathcal{S}$$ via $$re^{i\theta} \mapsto r \operatorname{Exp}(i\theta)$$. Here $$r > 0$$ and $$-\pi < \theta \le \pi$$. Then $$l$$ is holomorphic, and we define $$G = f \circ l$$. This called “lifting to the principal branch” throughout the SymPy documentation.

## Table Lookups and Inverse Mellin Transforms#

Suppose we are given an integrand $$f(x)$$ and are trying to rewrite it as a single G-function. To do this, we first split $$f(x)$$ into the form $$x^s g(x)$$ (where $$g(x)$$ is supposed to be simpler than $$f(x)$$). This is because multiplicative powers can be absorbed into the G-function later. This splitting is done by _split_mul(f, x). Then we assemble a tuple of functions that occur in $$f$$ (e.g. if $$f(x) = e^x \cos{x}$$, we would assemble the tuple $$(\cos, \exp)$$). This is done by the function _mytype(f, x). Next we index a lookup table (created using _create_lookup_table()) with this tuple. This (hopefully) yields a list of Meijer G-function formulae involving these functions, we then pattern-match all of them. If one fits, we were successful, otherwise not and we have to try something else.

Suppose now we want to rewrite as a product of two G-functions. To do this, we (try to) find all inequivalent ways of splitting $$f(x)$$ into a product $$f_1(x) f_2(x)$$. We could try these splittings in any order, but it is often a good idea to minimize (a) the number of powers occurring in $$f_i(x)$$ and (b) the number of different functions occurring in $$f_i(x)$$. Thus given e.g. $$f(x) = \sin{x}\, e^{x} \sin{2x}$$ we should try $$f_1(x) = \sin{x}\, \sin{2x}$$, $$f_2(x) = e^{x}$$ first. All of this is done by the function _mul_as_two_parts(f).

Finally, we can try a recursive Mellin transform technique. Since the Meijer G-function is defined essentially as a certain inverse mellin transform, if we want to write a function $$f(x)$$ as a G-function, we can compute its mellin transform $$F(s)$$. If $$F(s)$$ is in the right form, the G-function expression can be read off. This technique generalises many standard rewritings, e.g. $$e^{ax} e^{bx} = e^{(a + b) x}$$.

One twist is that some functions don’t have mellin transforms, even though they can be written as G-functions. This is true for example for $$f(x) = e^x \sin{x}$$ (the function grows too rapidly to have a mellin transform). However if the function is recognised to be analytic, then we can try to compute the mellin-transform of $$f(ax)$$ for a parameter $$a$$, and deduce the G-function expression by analytic continuation. (Checking for analyticity is easy. Since we can only deal with a certain subset of functions anyway, we only have to filter out those which are not analyitc.)

The function _rewrite_single does the table lookup and recursive mellin transform. The functions _rewrite1 and _rewrite2 respectively use above-mentioned helpers and _rewrite_single to rewrite their argument as respectively one or two G-functions.

## Applying the Integral Theorems#

If the integrand has been recast into G-functions, evaluating the integral is relatively easy. We first do some substitutions to reduce e.g. the exponent of the argument of the G-function to unity (see _rewrite_saxena_1 and _rewrite_saxena, respectively, for one or two G-functions). Next we go through a list of conditions under which the integral theorem applies. It can fail for basically two reasons: either the integral does not exist, or the manipulations in deriving the theorem may not be allowed (for more details, see this [BlogPost]).

Sometimes this can be remedied by reducing the argument of the G-functions involved. For example it is clear that the G-function representing $$e^z$$ is satisfies $$G(\operatorname{Exp}(2 \pi i)z) = G(z)$$ for all $$z \in \mathcal{S}$$. The function meijerg.get_period() can be used to discover this, and the function principal_branch(z, period) in functions/elementary/complexes.py can be used to exploit the information. This is done transparently by the integration code.

### The G-Function Integration Theorems#

This section intends to display in detail the definite integration theorems used in the code. The following two formulae go back to Meijer (In fact he proved more general formulae; indeed in the literature formulae are usually staded in more general form. However it is very easy to deduce the general formulae from the ones we give here. It seemed best to keep the theorems as simple as possible, since they are very complicated anyway.):

1. $\begin{split}\int_0^\infty G_{p, q}^{m, n} \left.\left(\begin{matrix} a_1, \cdots, a_p \\ b_1, \cdots, b_q \end{matrix} \right| \eta x \right) \mathrm{d}x = \frac{\prod_{j=1}^m \Gamma(b_j + 1) \prod_{j=1}^n \Gamma(-a_j)}{\eta \prod_{j=m+1}^q \Gamma(-b_j) \prod_{j=n+1}^p \Gamma(a_j + 1)}\end{split}$
2. $\begin{split}\int_0^\infty G_{u, v}^{s, t} \left.\left(\begin{matrix} c_1, \cdots, c_u \\ d_1, \cdots, d_v \end{matrix} \right| \sigma x \right) G_{p, q}^{m, n} \left.\left(\begin{matrix} a_1, \cdots, a_p \\ b_1, \cdots, b_q \end{matrix} \right| \omega x \right) \mathrm{d}x = G_{v+p, u+q}^{m+t, n+s} \left.\left( \begin{matrix} a_1, \cdots, a_n, -d_1, \cdots, -d_v, a_{n+1}, \cdots, a_p \\ b_1, \cdots, b_m, -c_1, \cdots, -c_u, b_{m+1}, \cdots, b_q \end{matrix} \right| \frac{\omega}{\sigma} \right)\end{split}$

The more interesting question is under what conditions these formulae are valid. Below we detail the conditions implemented in SymPy. They are an amalgamation of conditions found in [Prudnikov1990] and [Luke1969]; please let us know if you find any errors.

## Conditions of Convergence for Integral (1)#

We can without loss of generality assume $$p \le q$$, since the G-functions of indices $$m, n, p, q$$ and of indices $$n, m, q, p$$ can be related easily (see e.g. [Luke1969], section 5.3). We introduce the following notation:

$\begin{split}\xi = m + n - p \\ \delta = m + n - \frac{p + q}{2}\end{split}$
$\begin{split}C_3: -\Re(b_j) < 1 \text{ for } j=1, \ldots, m \\ 0 < -\Re(a_j) \text{ for } j=1, \ldots, n\end{split}$
$\begin{split}C_3^*: -\Re(b_j) < 1 \text{ for } j=1, \ldots, q \\ 0 < -\Re(a_j) \text{ for } j=1, \ldots, p\end{split}$
$C_4: -\Re(\delta) + \frac{q + 1 - p}{2} > q - p$

The convergence conditions will be detailed in several “cases”, numbered one to five. For later use it will be helpful to separate conditions “at infinity” from conditions “at zero”. By conditions “at infinity” we mean conditions that only depend on the behaviour of the integrand for large, positive values of $$x$$, whereas by conditions “at zero” we mean conditions that only depend on the behaviour of the integrand on $$(0, \epsilon)$$ for any $$\epsilon > 0$$. Since all our conditions are specified in terms of parameters of the G-functions, this distinction is not immediately visible. They are, however, of very distinct character mathematically; the conditions at infinity being in particular much harder to control.

In order for the integral theorem to be valid, conditions $$n$$ “at zero” and “at infinity” both have to be fulfilled, for some $$n$$.

These are the conditions “at infinity”:

1. $\delta > 0 \wedge |\arg(\eta)| < \delta \pi \wedge (A \vee B \vee C),$

where

$A = 1 \le n \wedge p < q \wedge 1 \le m$
$B = 1 \le p \wedge 1 \le m \wedge q = p+1 \wedge \neg (n = 0 \wedge m = p + 1 )$
$C = 1 \le n \wedge q = p \wedge |\arg(\eta)| \ne (\delta - 2k)\pi \text{ for } k = 0, 1, \ldots \left\lceil \frac{\delta}{2} \right\rceil.$
2. $n = 0 \wedge p + 1 \le m \wedge |\arg(\eta)| < \delta \pi$
3. $(p < q \wedge 1 \le m \wedge \delta > 0 \wedge |\arg(\eta)| = \delta \pi) \vee (p \le q - 2 \wedge \delta = 0 \wedge \arg(\eta) = 0)$
4. $p = q \wedge \delta = 0 \wedge \arg(\eta) = 0 \wedge \eta \ne 0 \wedge \Re\left(\sum_{j=1}^p b_j - a_j \right) < 0$
5. $\delta > 0 \wedge |\arg(\eta)| < \delta \pi$

And these are the conditions “at zero”:

1. $\eta \ne 0 \wedge C_3$
2. $C_3$
3. $C_3 \wedge C_4$
4. $C_3$
5. $C_3$

## Conditions of Convergence for Integral (2)#

We introduce the following notation:

$b^* = s + t - \frac{u + v}{2}$
$c^* = m + n - \frac{p + q}{2}$
$\rho = \sum_{j=1}^v d_j - \sum_{j=1}^u c_j + \frac{u - v}{2} + 1$
$\mu = \sum_{j=1}^q b_j - \sum_{j=1}^p a_j + \frac{p - q}{2} + 1$
$\phi = q - p - \frac{u - v}{2} + 1$
$\eta = 1 - (v - u) - \mu - \rho$
$\psi = \frac{\pi(q - m - n) + |\arg(\omega)|}{q - p}$
$\theta = \frac{\pi(v - s - t) + |\arg(\sigma)|)}{v - u}$
$\lambda_c = (q - p)|\omega|^{1/(q - p)} \cos{\psi} + (v - u)|\sigma|^{1/(v - u)} \cos{\theta}$
$\lambda_{s0}(c_1, c_2) = c_1 (q - p)|\omega|^{1/(q - p)} \sin{\psi} + c_2 (v - u)|\sigma|^{1/(v - u)} \sin{\theta}$
$\begin{split}\lambda_s = \begin{cases} \operatorname{\lambda_{s0}}\left(-1,-1\right) \operatorname{\lambda_{s0}}\left(1,1\right) & \text{for}\: \arg(\omega) = 0 \wedge \arg(\sigma) = 0 \\\operatorname{\lambda_{s0}}\left(\operatorname{sign}\left(\operatorname{\arg}\left(\omega\right)\right),-1\right) \operatorname{\lambda_{s0}}\left(\operatorname{sign}\left(\operatorname{\arg}\left(\omega\right)\right),1\right) & \text{for}\: \arg(\omega) \ne 0 \wedge \arg(\sigma) = 0 \\\operatorname{\lambda_{s0}}\left(-1,\operatorname{sign}\left(\operatorname{\arg}\left(\sigma\right)\right)\right) \operatorname{\lambda_{s0}}\left(1,\operatorname{sign}\left(\operatorname{\arg}\left(\sigma\right)\right)\right) & \text{for}\: \arg(\omega) = 0 \wedge \arg(\sigma) \ne 0) \\\operatorname{\lambda_{s0}}\left(\operatorname{sign}\left(\operatorname{\arg}\left(\omega\right)\right),\operatorname{sign}\left(\operatorname{\arg}\left(\sigma\right)\right)\right) & \text{otherwise} \end{cases}\end{split}$
$z_0 = \frac{\omega}{\sigma} e^{-i\pi (b^* + c^*)}$
$z_1 = \frac{\sigma}{\omega} e^{-i\pi (b^* + c^*)}$

The following conditions will be helpful:

$\begin{split}C_1: (a_i - b_j \notin \mathbb{Z}_{>0} \text{ for } i = 1, \ldots, n, j = 1, \ldots, m) \\ \wedge (c_i - d_j \notin \mathbb{Z}_{>0} \text{ for } i = 1, \ldots, t, j = 1, \ldots, s)\end{split}$
$C_2: \Re(1 + b_i + d_j) > 0 \text{ for } i = 1, \ldots, m, j = 1, \ldots, s$
$C_3: \Re(a_i + c_j) < 1 \text{ for } i = 1, \ldots, n, j = 1, \ldots, t$
$C_4: (p - q)\Re(c_i) - \Re(\mu) > -\frac{3}{2} \text{ for } i=1, \ldots, t$
$C_5: (p - q)\Re(1 + d_i) - \Re(\mu) > -\frac{3}{2} \text{ for } i=1, \ldots, s$
$C_6: (u - v)\Re(a_i) - \Re(\rho) > -\frac{3}{2} \text{ for } i=1, \ldots, n$
$C_7: (u - v)\Re(1 + b_i) - \Re(\rho) > -\frac{3}{2} \text{ for } i=1, \ldots, m$
$C_8: 0 < \lvert{\phi}\rvert + 2 \Re\left(\left(\mu -1\right) \left(- u + v\right) + \left(- p + q\right) \left(\rho -1\right) + \left(- p + q\right) \left(- u + v\right)\right)$
$C_9: 0 < \lvert{\phi}\rvert - 2 \Re\left(\left(\mu -1\right) \left(- u + v\right) + \left(- p + q\right) \left(\rho -1\right) + \left(- p + q\right) \left(- u + v\right)\right)$
$C_{10}: \lvert{\operatorname{arg}\left(\sigma\right)}\rvert < \pi b^{*}$
$C_{11}: \lvert{\operatorname{arg}\left(\sigma\right)}\rvert = \pi b^{*}$
$C_{12}: |\arg(\omega)| < c^*\pi$
$C_{13}: |\arg(\omega)| = c^*\pi$
$C_{14}^1: \left(z_0 \ne 1 \wedge |\arg(1 - z_0)| < \pi \right) \vee \left(z_0 = 1 \wedge \Re(\mu + \rho - u + v) < 1 \right)$
$C_{14}^2: \left(z_1 \ne 1 \wedge |\arg(1 - z_1)| < \pi \right) \vee \left(z_1 = 1 \wedge \Re(\mu + \rho - p + q) < 1 \right)$
$C_{14}: \phi = 0 \wedge b^* + c^* \le 1 \wedge (C_{14}^1 \vee C_{14}^2)$
$C_{15}: \lambda_c > 0 \vee (\lambda_c = 0 \wedge \lambda_s \ne 0 \wedge \Re(\eta) > -1) \vee (\lambda_c = 0 \wedge \lambda_s = 0 \wedge \Re(\eta) > 0)$
$C_{16}: \int_0^\infty G_{u, v}^{s, t}(\sigma x) \mathrm{d} x \text{ converges at infinity }$
$C_{17}: \int_0^\infty G_{p, q}^{m, n}(\omega x) \mathrm{d} x \text{ converges at infinity }$

Note that $$C_{16}$$ and $$C_{17}$$ are the reason we split the convergence conditions for integral (1).

With this notation established, the implemented convergence conditions can be enumerated as follows:

1. $m n s t \neq 0 \wedge 0 < b^{*} \wedge 0 < c^{*} \wedge C_{1} \wedge C_{2} \wedge C_{3} \wedge C_{10} \wedge C_{12}$
2. $u = v \wedge b^{*} = 0 \wedge 0 < c^{*} \wedge 0 < \sigma \wedge \Re{\rho} < 1 \wedge C_{1} \wedge C_{2} \wedge C_{3} \wedge C_{12}$
3. $p = q \wedge u = v \wedge b^{*} = 0 \wedge c^{*} = 0 \wedge 0 < \sigma \wedge 0 < \omega \wedge \Re{\mu} < 1 \wedge \Re{\rho} < 1 \wedge \sigma \neq \omega \wedge C_{1} \wedge C_{2} \wedge C_{3}$
4. $p = q \wedge u = v \wedge b^{*} = 0 \wedge c^{*} = 0 \wedge 0 < \sigma \wedge 0 < \omega \wedge \Re\left(\mu + \rho\right) < 1 \wedge \omega \neq \sigma \wedge C_{1} \wedge C_{2} \wedge C_{3}$
5. $p = q \wedge u = v \wedge b^{*} = 0 \wedge c^{*} = 0 \wedge 0 < \sigma \wedge 0 < \omega \wedge \Re\left(\mu + \rho\right) < 1 \wedge \omega \neq \sigma \wedge C_{1} \wedge C_{2} \wedge C_{3}$
6. $q < p \wedge 0 < s \wedge 0 < b^{*} \wedge 0 \leq c^{*} \wedge C_{1} \wedge C_{2} \wedge C_{3} \wedge C_{5} \wedge C_{10} \wedge C_{13}$
7. $p < q \wedge 0 < t \wedge 0 < b^{*} \wedge 0 \leq c^{*} \wedge C_{1} \wedge C_{2} \wedge C_{3} \wedge C_{4} \wedge C_{10} \wedge C_{13}$
8. $v < u \wedge 0 < m \wedge 0 < c^{*} \wedge 0 \leq b^{*} \wedge C_{1} \wedge C_{2} \wedge C_{3} \wedge C_{7} \wedge C_{11} \wedge C_{12}$
9. $u < v \wedge 0 < n \wedge 0 < c^{*} \wedge 0 \leq b^{*} \wedge C_{1} \wedge C_{2} \wedge C_{3} \wedge C_{6} \wedge C_{11} \wedge C_{12}$
10. $q < p \wedge u = v \wedge b^{*} = 0 \wedge 0 \leq c^{*} \wedge 0 < \sigma \wedge \Re{\rho} < 1 \wedge C_{1} \wedge C_{2} \wedge C_{3} \wedge C_{5} \wedge C_{13}$
11. $p < q \wedge u = v \wedge b^{*} = 0 \wedge 0 \leq c^{*} \wedge 0 < \sigma \wedge \Re{\rho} < 1 \wedge C_{1} \wedge C_{2} \wedge C_{3} \wedge C_{4} \wedge C_{13}$
12. $p = q \wedge v < u \wedge 0 \leq b^{*} \wedge c^{*} = 0 \wedge 0 < \omega \wedge \Re{\mu} < 1 \wedge C_{1} \wedge C_{2} \wedge C_{3} \wedge C_{7} \wedge C_{11}$
13. $p = q \wedge u < v \wedge 0 \leq b^{*} \wedge c^{*} = 0 \wedge 0 < \omega \wedge \Re{\mu} < 1 \wedge C_{1} \wedge C_{2} \wedge C_{3} \wedge C_{6} \wedge C_{11}$
14. $p < q \wedge v < u \wedge 0 \leq b^{*} \wedge 0 \leq c^{*} \wedge C_{1} \wedge C_{2} \wedge C_{3} \wedge C_{4} \wedge C_{7} \wedge C_{11} \wedge C_{13}$
15. $q < p \wedge u < v \wedge 0 \leq b^{*} \wedge 0 \leq c^{*} \wedge C_{1} \wedge C_{2} \wedge C_{3} \wedge C_{5} \wedge C_{6} \wedge C_{11} \wedge C_{13}$
16. $q < p \wedge v < u \wedge 0 \leq b^{*} \wedge 0 \leq c^{*} \wedge C_{1} \wedge C_{2} \wedge C_{3} \wedge C_{5} \wedge C_{7} \wedge C_{8} \wedge C_{11} \wedge C_{13} \wedge C_{14}$
17. $p < q \wedge u < v \wedge 0 \leq b^{*} \wedge 0 \leq c^{*} \wedge C_{1} \wedge C_{2} \wedge C_{3} \wedge C_{4} \wedge C_{6} \wedge C_{9} \wedge C_{11} \wedge C_{13} \wedge C_{14}$
18. $t = 0 \wedge 0 < s \wedge 0 < b^{*} \wedge 0 < \phi \wedge C_{1} \wedge C_{2} \wedge C_{10}$
19. $s = 0 \wedge 0 < t \wedge 0 < b^{*} \wedge \phi < 0 \wedge C_{1} \wedge C_{3} \wedge C_{10}$
20. $n = 0 \wedge 0 < m \wedge 0 < c^{*} \wedge \phi < 0 \wedge C_{1} \wedge C_{2} \wedge C_{12}$
21. $m = 0 \wedge 0 < n \wedge 0 < c^{*} \wedge 0 < \phi \wedge C_{1} \wedge C_{3} \wedge C_{12}$
22. $s t = 0 \wedge 0 < b^{*} \wedge 0 < c^{*} \wedge C_{1} \wedge C_{2} \wedge C_{3} \wedge C_{10} \wedge C_{12}$
23. $m n = 0 \wedge 0 < b^{*} \wedge 0 < c^{*} \wedge C_{1} \wedge C_{2} \wedge C_{3} \wedge C_{10} \wedge C_{12}$
24. $p < m + n \wedge t = 0 \wedge \phi = 0 \wedge 0 < s \wedge 0 < b^{*} \wedge c^{*} < 0 \wedge \lvert{\operatorname{arg}\left(\omega\right)}\rvert < \pi \left(m + n - p + 1\right) \wedge C_{1} \wedge C_{2} \wedge C_{10} \wedge C_{14} \wedge C_{15}$
25. $q < m + n \wedge s = 0 \wedge \phi = 0 \wedge 0 < t \wedge 0 < b^{*} \wedge c^{*} < 0 \wedge \lvert{\operatorname{arg}\left(\omega\right)}\rvert < \pi \left(m + n - q + 1\right) \wedge C_{1} \wedge C_{3} \wedge C_{10} \wedge C_{14} \wedge C_{15}$
26. $p = q -1 \wedge t = 0 \wedge \phi = 0 \wedge 0 < s \wedge 0 < b^{*} \wedge 0 \leq c^{*} \wedge \pi c^{*} < \lvert{\operatorname{arg}\left(\omega\right)}\rvert \wedge C_{1} \wedge C_{2} \wedge C_{10} \wedge C_{14} \wedge C_{15}$
27. $p = q + 1 \wedge s = 0 \wedge \phi = 0 \wedge 0 < t \wedge 0 < b^{*} \wedge 0 \leq c^{*} \wedge \pi c^{*} < \lvert{\operatorname{arg}\left(\omega\right)}\rvert \wedge C_{1} \wedge C_{3} \wedge C_{10} \wedge C_{14} \wedge C_{15}$
28. $p < q -1 \wedge t = 0 \wedge \phi = 0 \wedge 0 < s \wedge 0 < b^{*} \wedge 0 \leq c^{*} \wedge \pi c^{*} < \lvert{\operatorname{arg}\left(\omega\right)}\rvert \wedge \lvert{\operatorname{arg}\left(\omega\right)}\rvert < \pi \left(m + n - p + 1\right) \wedge C_{1} \wedge C_{2} \wedge C_{10} \wedge C_{14} \wedge C_{15}$
29. $q + 1 < p \wedge s = 0 \wedge \phi = 0 \wedge 0 < t \wedge 0 < b^{*} \wedge 0 \leq c^{*} \wedge \pi c^{*} < \lvert{\operatorname{arg}\left(\omega\right)}\rvert \wedge \lvert{\operatorname{arg}\left(\omega\right)}\rvert < \pi \left(m + n - q + 1 \right) \wedge C_{1} \wedge C_{3} \wedge C_{10} \wedge C_{14} \wedge C_{15}$
30. $n = 0 \wedge \phi = 0 \wedge 0 < s + t \wedge 0 < m \wedge 0 < c^{*} \wedge b^{*} < 0 \wedge \lvert{\operatorname{arg}\left(\sigma\right)}\rvert < \pi \left(s + t - u + 1\right) \wedge C_{1} \wedge C_{2} \wedge C_{12} \wedge C_{14} \wedge C_{15}$
31. $m = 0 \wedge \phi = 0 \wedge v < s + t \wedge 0 < n \wedge 0 < c^{*} \wedge b^{*} < 0 \wedge \lvert{\operatorname{arg}\left(\sigma\right)}\rvert < \pi \left(s + t - v + 1\right) \wedge C_{1} \wedge C_{3} \wedge C_{12} \wedge C_{14} \wedge C_{15}$
32. $n = 0 \wedge \phi = 0 \wedge u = v -1 \wedge 0 < m \wedge 0 < c^{*} \wedge 0 \leq b^{*} \wedge \pi b^{*} < \lvert{\operatorname{arg}\left(\sigma\right)}\rvert \wedge \lvert{\operatorname{arg}\left(\sigma\right)}\rvert < \pi \left(b^{*} + 1\right) \wedge C_{1} \wedge C_{2} \wedge C_{12} \wedge C_{14} \wedge C_{15}$
33. $m = 0 \wedge \phi = 0 \wedge u = v + 1 \wedge 0 < n \wedge 0 < c^{*} \wedge 0 \leq b^{*} \wedge \pi b^{*} < \lvert{\operatorname{arg}\left(\sigma\right)}\rvert \wedge \lvert{\operatorname{arg}\left(\sigma\right)}\rvert < \pi \left(b^{*} + 1\right) \wedge C_{1} \wedge C_{3} \wedge C_{12} \wedge C_{14} \wedge C_{15}$
34. $n = 0 \wedge \phi = 0 \wedge u < v -1 \wedge 0 < m \wedge 0 < c^{*} \wedge 0 \leq b^{*} \wedge \pi b^{*} < \lvert{\operatorname{arg}\left(\sigma\right)}\rvert \wedge \lvert{\operatorname{arg}\left(\sigma\right)}\rvert < \pi \left(s + t - u + 1\right) \wedge C_{1} \wedge C_{2} \wedge C_{12} \wedge C_{14} \wedge C_{15}$
35. $m = 0 \wedge \phi = 0 \wedge v + 1 < u \wedge 0 < n \wedge 0 < c^{*} \wedge 0 \leq b^{*} \wedge \pi b^{*} < \lvert{\operatorname{arg}\left(\sigma\right)}\rvert \wedge \lvert{\operatorname{arg}\left(\sigma\right)}\rvert < \pi \left(s + t - v + 1 \right) \wedge C_{1} \wedge C_{3} \wedge C_{12} \wedge C_{14} \wedge C_{15}$
36. $C_{17} \wedge t = 0 \wedge u < s \wedge 0 < b^{*} \wedge C_{10} \wedge C_{1} \wedge C_{2} \wedge C_{3}$
37. $C_{17} \wedge s = 0 \wedge v < t \wedge 0 < b^{*} \wedge C_{10} \wedge C_{1} \wedge C_{2} \wedge C_{3}$
38. $C_{16} \wedge n = 0 \wedge p < m \wedge 0 < c^{*} \wedge C_{12} \wedge C_{1} \wedge C_{2} \wedge C_{3}$
39. $C_{16} \wedge m = 0 \wedge q < n \wedge 0 < c^{*} \wedge C_{12} \wedge C_{1} \wedge C_{2} \wedge C_{3}$

### The Inverse Laplace Transform of a G-function#

The inverse laplace transform of a Meijer G-function can be expressed as another G-function. This is a fairly versatile method for computing this transform. However, I could not find the details in the literature, so I work them out here. In [Luke1969], section 5.6.3, there is a formula for the inverse Laplace transform of a G-function of argument $$bz$$, and convergence conditions are also given. However, we need a formula for argument $$bz^a$$ for rational $$a$$.

We are asked to compute

$f(t) = \frac{1}{2\pi i} \int_{c-i\infty}^{c+i\infty} e^{zt} G(bz^a) \mathrm{d}z,$

for positive real $$t$$. Three questions arise:

1. When does this integral converge?

2. How can we compute the integral?

3. When is our computation valid?

## How to compute the integral#

We shall work formally for now. Denote by $$\Delta(s)$$ the product of gamma functions appearing in the definition of $$G$$, so that

$G(z) = \frac{1}{2\pi i} \int_L \Delta(s) z^s \mathrm{d}s.$

Thus

$f(t) = \frac{1}{(2\pi i)^2} \int_{c - i\infty}^{c + i\infty} \int_L e^{zt} \Delta(s) b^s z^{as} \mathrm{d}s \mathrm{d}z.$

We interchange the order of integration to get

$f(t) = \frac{1}{2\pi i} \int_L b^s \Delta(s) \int_{c-i\infty}^{c+i\infty} e^{zt} z^{as} \frac{\mathrm{d}z}{2\pi i} \mathrm{d}s.$

The inner integral is easily seen to be $$\frac{1}{\Gamma(-as)} \frac{1}{t^{1+as}}$$. (Using Cauchy’s theorem and Jordan’s lemma deform the contour to run from $$-\infty$$ to $$-\infty$$, encircling $$0$$ once in the negative sense. For $$as$$ real and greater than one, this contour can be pushed onto the negative real axis and the integral is recognised as a product of a sine and a gamma function. The formula is then proved using the functional equation of the gamma function, and extended to the entire domain of convergence of the original integral by appealing to analytic continuation.) Hence we find

$f(t) = \frac{1}{t} \frac{1}{2\pi i} \int_L \Delta(s) \frac{1}{\Gamma(-as)} \left(\frac{b}{t^a}\right)^s \mathrm{d}s,$

which is a so-called Fox H function (of argument $$\frac{b}{t^a}$$). For rational $$a$$, this can be expressed as a Meijer G-function using the gamma function multiplication theorem.

## When this computation is valid#

There are a number of obstacles in this computation. Interchange of integrals is only valid if all integrals involved are absolutely convergent. In particular the inner integral has to converge. Also, for our identification of the final integral as a Fox H / Meijer G-function to be correct, the poles of the newly obtained gamma function must be separated properly.

It is easy to check that the inner integral converges absolutely for $$\Re(as) < -1$$. Thus the contour $$L$$ has to run left of the line $$\Re(as) = -1$$. Under this condition, the poles of the newly-introduced gamma function are separated properly.

It remains to observe that the Meijer G-function is an analytic, unbranched function of its parameters, and of the coefficient $$b$$. Hence so is $$f(t)$$. Thus the final computation remains valid as long as the initial integral converges, and if there exists a changed set of parameters where the computation is valid. If we assume w.l.o.g. that $$a > 0$$, then the latter condition is fulfilled if $$G$$ converges along contours (2) or (3) of [Luke1969], section 5.2, i.e. either $$\delta \ge \frac{a}{2}$$ or $$p \ge 1, p \ge q$$.

## When the integral exists#

Using [Luke1969], section 5.10, for any given meijer G-function we can find a dominant term of the form $$z^a e^{bz^c}$$ (although this expression might not be the best possible, because of cancellation).

We must thus investigate

$\lim_{T \to \infty} \int_{c-iT}^{c+iT} e^{zt} z^a e^{bz^c} \mathrm{d}z.$

(This principal value integral is the exact statement used in the Laplace inversion theorem.) We write $$z = c + i \tau$$. Then $$arg(z) \to \pm \frac{\pi}{2}$$, and so $$e^{zt} \sim e^{it \tau}$$ (where $$\sim$$ shall always mean “asymptotically equivalent up to a positive real multiplicative constant”). Also $$z^{x + iy} \sim |\tau|^x e^{i y \log{|\tau|}} e^{\pm x i \frac{\pi}{2}}.$$

Set $$\omega_{\pm} = b e^{\pm i \Re(c) \frac{\pi}{2}}$$. We have three cases:

1. $$b=0$$ or $$\Re(c) \le 0$$. In this case the integral converges if $$\Re(a) \le -1$$.

2. $$b \ne 0$$, $$\Im(c) = 0$$, $$\Re(c) > 0$$. In this case the integral converges if $$\Re(\omega_{\pm}) < 0$$.

3. $$b \ne 0$$, $$\Im(c) = 0$$, $$\Re(c) > 0$$, $$\Re(\omega_{\pm}) \le 0$$, and at least one of $$\Re(\omega_{\pm}) = 0$$. Here the same condition as in (1) applies.

### Implemented G-Function Formulae#

An important part of the algorithm is a table expressing various functions as Meijer G-functions. This is essentially a table of Mellin Transforms in disguise. The following automatically generated table shows the formulae currently implemented in SymPy. An entry “generated” means that the corresponding G-function has a variable number of parameters. This table is intended to shrink in future, when the algorithm’s capabilities of deriving new formulae improve. Of course it has to grow whenever a new class of special functions is to be dealt with.

Elementary functions:

$\begin{split}a = a {G_{1, 1}^{1, 0}\left(\begin{matrix} & 1 \\0 & \end{matrix} \middle| {z} \right)} + a {G_{1, 1}^{0, 1}\left(\begin{matrix} 1 & \\ & 0 \end{matrix} \middle| {z} \right)}\end{split}$
$\begin{split}\left(z^{q} p + b\right)^{- a} = \frac{b^{- a} {G_{1, 1}^{1, 1}\left(\begin{matrix} 1 - a & \\0 & \end{matrix} \middle| {\frac{z^{q} p}{b}} \right)}}{\Gamma\left(a\right)}\end{split}$
$\begin{split}\frac{- b^{a} + \left(z^{q} p\right)^{a}}{z^{q} p - b} = \frac{b^{a - 1} {G_{2, 2}^{2, 2}\left(\begin{matrix} 0, a & \\0, a & \end{matrix} \middle| {\frac{z^{q} p}{b}} \right)} \sin{\left(\pi a \right)}}{\pi}\end{split}$
$\begin{split}\left(a + \sqrt{z^{q} p + a^{2}}\right)^{b} = - \frac{a^{b} b {G_{2, 2}^{1, 2}\left(\begin{matrix} \frac{b}{2} + \frac{1}{2}, \frac{b}{2} + 1 & \\0 & b \end{matrix} \middle| {\frac{z^{q} p}{a^{2}}} \right)}}{2 \sqrt{\pi}}\end{split}$
$\begin{split}\left(- a + \sqrt{z^{q} p + a^{2}}\right)^{b} = \frac{a^{b} b {G_{2, 2}^{1, 2}\left(\begin{matrix} \frac{b}{2} + \frac{1}{2}, \frac{b}{2} + 1 & \\b & 0 \end{matrix} \middle| {\frac{z^{q} p}{a^{2}}} \right)}}{2 \sqrt{\pi}}\end{split}$
$\begin{split}\frac{\left(a + \sqrt{z^{q} p + a^{2}}\right)^{b}}{\sqrt{z^{q} p + a^{2}}} = \frac{a^{b - 1} {G_{2, 2}^{1, 2}\left(\begin{matrix} \frac{b}{2} + \frac{1}{2}, \frac{b}{2} & \\0 & b \end{matrix} \middle| {\frac{z^{q} p}{a^{2}}} \right)}}{\sqrt{\pi}}\end{split}$
$\begin{split}\frac{\left(- a + \sqrt{z^{q} p + a^{2}}\right)^{b}}{\sqrt{z^{q} p + a^{2}}} = \frac{a^{b - 1} {G_{2, 2}^{1, 2}\left(\begin{matrix} \frac{b}{2} + \frac{1}{2}, \frac{b}{2} & \\b & 0 \end{matrix} \middle| {\frac{z^{q} p}{a^{2}}} \right)}}{\sqrt{\pi}}\end{split}$
$\begin{split}\left(z^{\frac{q}{2}} \sqrt{p} + \sqrt{z^{q} p + a}\right)^{b} = - \frac{a^{\frac{b}{2}} b {G_{2, 2}^{2, 1}\left(\begin{matrix} \frac{b}{2} + 1 & 1 - \frac{b}{2} \\0, \frac{1}{2} & \end{matrix} \middle| {\frac{z^{q} p}{a}} \right)}}{2 \sqrt{\pi}}\end{split}$
$\begin{split}\left(- z^{\frac{q}{2}} \sqrt{p} + \sqrt{z^{q} p + a}\right)^{b} = \frac{a^{\frac{b}{2}} b {G_{2, 2}^{2, 1}\left(\begin{matrix} 1 - \frac{b}{2} & \frac{b}{2} + 1 \\0, \frac{1}{2} & \end{matrix} \middle| {\frac{z^{q} p}{a}} \right)}}{2 \sqrt{\pi}}\end{split}$
$\begin{split}\frac{\left(z^{\frac{q}{2}} \sqrt{p} + \sqrt{z^{q} p + a}\right)^{b}}{\sqrt{z^{q} p + a}} = \frac{a^{\frac{b}{2} - \frac{1}{2}} {G_{2, 2}^{2, 1}\left(\begin{matrix} \frac{b}{2} + \frac{1}{2} & \frac{1}{2} - \frac{b}{2} \\0, \frac{1}{2} & \end{matrix} \middle| {\frac{z^{q} p}{a}} \right)}}{\sqrt{\pi}}\end{split}$
$\begin{split}\frac{\left(- z^{\frac{q}{2}} \sqrt{p} + \sqrt{z^{q} p + a}\right)^{b}}{\sqrt{z^{q} p + a}} = \frac{a^{\frac{b}{2} - \frac{1}{2}} {G_{2, 2}^{2, 1}\left(\begin{matrix} \frac{1}{2} - \frac{b}{2} & \frac{b}{2} + \frac{1}{2} \\0, \frac{1}{2} & \end{matrix} \middle| {\frac{z^{q} p}{a}} \right)}}{\sqrt{\pi}}\end{split}$

Functions involving $$\left|{z^{q} p - b}\right|$$:

$\begin{split}\left|{z^{q} p - b}\right|^{- a} = 2 {G_{2, 2}^{1, 1}\left(\begin{matrix} 1 - a & \frac{1}{2} - \frac{a}{2} \\0 & \frac{1}{2} - \frac{a}{2} \end{matrix} \middle| {\frac{z^{q} p}{b}} \right)} \sin{\left(\frac{\pi a}{2} \right)} \left|{b}\right|^{- a} \Gamma\left(1 - a\right),\text{ if } \operatorname{re}{\left(a\right)} < 1\end{split}$

Functions involving $$\operatorname{Chi}\left(z^{q} p\right)$$:

$\begin{split}\operatorname{Chi}\left(z^{q} p\right) = - \frac{\pi^{\frac{3}{2}} {G_{2, 4}^{2, 0}\left(\begin{matrix} & \frac{1}{2}, 1 \\0, 0 & \frac{1}{2}, \frac{1}{2} \end{matrix} \middle| {\frac{z^{2 q} p^{2}}{4}} \right)}}{2}\end{split}$

Functions involving $$\operatorname{Ci}{\left(z^{q} p \right)}$$:

$\begin{split}\operatorname{Ci}{\left(z^{q} p \right)} = - \frac{\sqrt{\pi} {G_{1, 3}^{2, 0}\left(\begin{matrix} & 1 \\0, 0 & \frac{1}{2} \end{matrix} \middle| {\frac{z^{2 q} p^{2}}{4}} \right)}}{2}\end{split}$

Functions involving $$\operatorname{Ei}{\left(z^{q} p \right)}$$:

$\begin{split}\operatorname{Ei}{\left(z^{q} p \right)} = - i \pi {G_{1, 1}^{1, 0}\left(\begin{matrix} & 1 \\0 & \end{matrix} \middle| {z} \right)} - {G_{1, 2}^{2, 0}\left(\begin{matrix} & 1 \\0, 0 & \end{matrix} \middle| {z^{q} p e^{i \pi}} \right)} - i \pi {G_{1, 1}^{0, 1}\left(\begin{matrix} 1 & \\ & 0 \end{matrix} \middle| {z} \right)}\end{split}$

Functions involving $$\theta\left(z^{q} p - b\right)$$:

$\begin{split}\left(z^{q} p - b\right)^{a - 1} \theta\left(z^{q} p - b\right) = b^{a - 1} {G_{1, 1}^{0, 1}\left(\begin{matrix} a & \\ & 0 \end{matrix} \middle| {\frac{z^{q} p}{b}} \right)} \Gamma\left(a\right),\text{ if } b > 0\end{split}$
$\begin{split}\left(- z^{q} p + b\right)^{a - 1} \theta\left(- z^{q} p + b\right) = b^{a - 1} {G_{1, 1}^{1, 0}\left(\begin{matrix} & a \\0 & \end{matrix} \middle| {\frac{z^{q} p}{b}} \right)} \Gamma\left(a\right),\text{ if } b > 0\end{split}$
$\begin{split}\left(z^{q} p - b\right)^{a - 1} \theta\left(z - \left(\frac{b}{p}\right)^{\frac{1}{q}}\right) = b^{a - 1} {G_{1, 1}^{0, 1}\left(\begin{matrix} a & \\ & 0 \end{matrix} \middle| {\frac{z^{q} p}{b}} \right)} \Gamma\left(a\right),\text{ if } b > 0\end{split}$
$\begin{split}\left(- z^{q} p + b\right)^{a - 1} \theta\left(- z + \left(\frac{b}{p}\right)^{\frac{1}{q}}\right) = b^{a - 1} {G_{1, 1}^{1, 0}\left(\begin{matrix} & a \\0 & \end{matrix} \middle| {\frac{z^{q} p}{b}} \right)} \Gamma\left(a\right),\text{ if } b > 0\end{split}$

Functions involving $$\theta\left(- z^{q} p + 1\right)$$, $$\log{\left(z^{q} p \right)}$$:

$\log{\left(z^{q} p \right)}^{n} \theta\left(- z^{q} p + 1\right) = \text{generated}$
$\log{\left(z^{q} p \right)}^{n} \theta\left(z^{q} p - 1\right) = \text{generated}$

Functions involving $$\operatorname{Shi}{\left(z^{q} p \right)}$$:

$\begin{split}\operatorname{Shi}{\left(z^{q} p \right)} = \frac{z^{q} \sqrt{\pi} p {G_{1, 3}^{1, 1}\left(\begin{matrix} \frac{1}{2} & \\0 & - \frac{1}{2}, - \frac{1}{2} \end{matrix} \middle| {\frac{z^{2 q} p^{2} e^{i \pi}}{4}} \right)}}{4}\end{split}$

Functions involving $$\operatorname{Si}{\left(z^{q} p \right)}$$:

$\begin{split}\operatorname{Si}{\left(z^{q} p \right)} = \frac{\sqrt{\pi} {G_{1, 3}^{1, 1}\left(\begin{matrix} 1 & \\\frac{1}{2} & 0, 0 \end{matrix} \middle| {\frac{z^{2 q} p^{2}}{4}} \right)}}{2}\end{split}$

Functions involving $$I_{a}\left(z^{q} p\right)$$:

$\begin{split}I_{a}\left(z^{q} p\right) = \pi {G_{1, 3}^{1, 0}\left(\begin{matrix} & \frac{a}{2} + \frac{1}{2} \\\frac{a}{2} & - \frac{a}{2}, \frac{a}{2} + \frac{1}{2} \end{matrix} \middle| {\frac{z^{2 q} p^{2}}{4}} \right)}\end{split}$

Functions involving $$J_{a}\left(z^{q} p\right)$$:

$\begin{split}J_{a}\left(z^{q} p\right) = {G_{0, 2}^{1, 0}\left(\begin{matrix} & \\\frac{a}{2} & - \frac{a}{2} \end{matrix} \middle| {\frac{z^{2 q} p^{2}}{4}} \right)}\end{split}$

Functions involving $$K_{a}\left(z^{q} p\right)$$:

$\begin{split}K_{a}\left(z^{q} p\right) = \frac{{G_{0, 2}^{2, 0}\left(\begin{matrix} & \\\frac{a}{2}, - \frac{a}{2} & \end{matrix} \middle| {\frac{z^{2 q} p^{2}}{4}} \right)}}{2}\end{split}$

Functions involving $$Y_{a}\left(z^{q} p\right)$$:

$\begin{split}Y_{a}\left(z^{q} p\right) = {G_{1, 3}^{2, 0}\left(\begin{matrix} & - \frac{a}{2} - \frac{1}{2} \\\frac{a}{2}, - \frac{a}{2} & - \frac{a}{2} - \frac{1}{2} \end{matrix} \middle| {\frac{z^{2 q} p^{2}}{4}} \right)}\end{split}$

Functions involving $$\cos{\left(z^{q} p \right)}$$:

$\begin{split}\cos{\left(z^{q} p \right)} = \sqrt{\pi} {G_{0, 2}^{1, 0}\left(\begin{matrix} & \\0 & \frac{1}{2} \end{matrix} \middle| {\frac{z^{2 q} p^{2}}{4}} \right)}\end{split}$

Functions involving $$\cosh{\left(z^{q} p \right)}$$:

$\begin{split}\cosh{\left(z^{q} p \right)} = \pi^{\frac{3}{2}} {G_{1, 3}^{1, 0}\left(\begin{matrix} & \frac{1}{2} \\0 & \frac{1}{2}, \frac{1}{2} \end{matrix} \middle| {\frac{z^{2 q} p^{2}}{4}} \right)}\end{split}$

Functions involving $$E\left(z^{q} p\right)$$:

$\begin{split}E\left(z^{q} p\right) = - \frac{{G_{2, 2}^{1, 2}\left(\begin{matrix} \frac{1}{2}, \frac{3}{2} & \\0 & 0 \end{matrix} \middle| {- z^{q} p} \right)}}{4}\end{split}$

Functions involving $$K\left(z^{q} p\right)$$:

$\begin{split}K\left(z^{q} p\right) = \frac{{G_{2, 2}^{1, 2}\left(\begin{matrix} \frac{1}{2}, \frac{1}{2} & \\0 & 0 \end{matrix} \middle| {- z^{q} p} \right)}}{2}\end{split}$

Functions involving $$\operatorname{erf}{\left(z^{q} p \right)}$$:

$\begin{split}\operatorname{erf}{\left(z^{q} p \right)} = \frac{{G_{1, 2}^{1, 1}\left(\begin{matrix} 1 & \\\frac{1}{2} & 0 \end{matrix} \middle| {z^{2 q} p^{2}} \right)}}{\sqrt{\pi}}\end{split}$

Functions involving $$\operatorname{erfc}{\left(z^{q} p \right)}$$:

$\begin{split}\operatorname{erfc}{\left(z^{q} p \right)} = \frac{{G_{1, 2}^{2, 0}\left(\begin{matrix} & 1 \\0, \frac{1}{2} & \end{matrix} \middle| {z^{2 q} p^{2}} \right)}}{\sqrt{\pi}}\end{split}$

Functions involving $$\operatorname{erfi}{\left(z^{q} p \right)}$$:

$\begin{split}\operatorname{erfi}{\left(z^{q} p \right)} = \frac{z^{q} p {G_{1, 2}^{1, 1}\left(\begin{matrix} \frac{1}{2} & \\0 & - \frac{1}{2} \end{matrix} \middle| {- z^{2 q} p^{2}} \right)}}{\sqrt{\pi}}\end{split}$

Functions involving $$e^{z^{q} p e^{i \pi}}$$:

$\begin{split}e^{z^{q} p e^{i \pi}} = {G_{0, 1}^{1, 0}\left(\begin{matrix} & \\0 & \end{matrix} \middle| {z^{q} p} \right)}\end{split}$

Functions involving $$\operatorname{E}_{a}\left(z^{q} p\right)$$:

$\begin{split}\operatorname{E}_{a}\left(z^{q} p\right) = {G_{1, 2}^{2, 0}\left(\begin{matrix} & a \\a - 1, 0 & \end{matrix} \middle| {z^{q} p} \right)}\end{split}$

Functions involving $$C\left(z^{q} p\right)$$:

$\begin{split}C\left(z^{q} p\right) = \frac{{G_{1, 3}^{1, 1}\left(\begin{matrix} 1 & \\\frac{1}{4} & 0, \frac{3}{4} \end{matrix} \middle| {\frac{z^{4 q} \pi^{2} p^{4}}{16}} \right)}}{2}\end{split}$

Functions involving $$S\left(z^{q} p\right)$$:

$\begin{split}S\left(z^{q} p\right) = \frac{{G_{1, 3}^{1, 1}\left(\begin{matrix} 1 & \\\frac{3}{4} & 0, \frac{1}{4} \end{matrix} \middle| {\frac{z^{4 q} \pi^{2} p^{4}}{16}} \right)}}{2}\end{split}$

Functions involving $$\log{\left(z^{q} p \right)}$$:

$\log{\left(z^{q} p \right)}^{n} = \text{generated}$
$\begin{split}\log{\left(z^{q} p + a \right)} = {G_{1, 1}^{1, 0}\left(\begin{matrix} & 1 \\0 & \end{matrix} \middle| {z} \right)} \log{\left(a \right)} + {G_{1, 1}^{0, 1}\left(\begin{matrix} 1 & \\ & 0 \end{matrix} \middle| {z} \right)} \log{\left(a \right)} + {G_{2, 2}^{1, 2}\left(\begin{matrix} 1, 1 & \\1 & 0 \end{matrix} \middle| {\frac{z^{q} p}{a}} \right)}\end{split}$
$\begin{split}\log{\left(\left|{z^{q} p - a}\right| \right)} = {G_{1, 1}^{1, 0}\left(\begin{matrix} & 1 \\0 & \end{matrix} \middle| {z} \right)} \log{\left(\left|{a}\right| \right)} + {G_{1, 1}^{0, 1}\left(\begin{matrix} 1 & \\ & 0 \end{matrix} \middle| {z} \right)} \log{\left(\left|{a}\right| \right)} + \pi {G_{3, 3}^{1, 2}\left(\begin{matrix} 1, 1 & \frac{1}{2} \\1 & 0, \frac{1}{2} \end{matrix} \middle| {\frac{z^{q} p}{a}} \right)}\end{split}$

Functions involving $$\sin{\left(z^{q} p \right)}$$:

$\begin{split}\sin{\left(z^{q} p \right)} = \sqrt{\pi} {G_{0, 2}^{1, 0}\left(\begin{matrix} & \\\frac{1}{2} & 0 \end{matrix} \middle| {\frac{z^{2 q} p^{2}}{4}} \right)}\end{split}$

Functions involving $$\operatorname{sinc}{\left(z^{q} p \right)}$$:

$\begin{split}\operatorname{sinc}{\left(z^{q} p \right)} = \frac{\sqrt{\pi} {G_{0, 2}^{1, 0}\left(\begin{matrix} & \\0 & - \frac{1}{2} \end{matrix} \middle| {\frac{z^{2 q} p^{2}}{4}} \right)}}{2}\end{split}$

Functions involving $$\sinh{\left(z^{q} p \right)}$$:

$\begin{split}\sinh{\left(z^{q} p \right)} = \pi^{\frac{3}{2}} {G_{1, 3}^{1, 0}\left(\begin{matrix} & 1 \\\frac{1}{2} & 1, 0 \end{matrix} \middle| {\frac{z^{2 q} p^{2}}{4}} \right)}\end{split}$

### Internal API Reference#

Integrate functions by rewriting them as Meijer G-functions.

There are three user-visible functions that can be used by other parts of the sympy library to solve various integration problems:

• meijerint_indefinite

• meijerint_definite

• meijerint_inversion

They can be used to compute, respectively, indefinite integrals, definite integrals over intervals of the real line, and inverse laplace-type integrals (from c-I*oo to c+I*oo). See the respective docstrings for details.

The main references for this are:

[L] Luke, Y. L. (1969), The Special Functions and Their Approximations,

Volume 1

[R] Kelly B. Roach. Meijer G Function Representations.

In: Proceedings of the 1997 International Symposium on Symbolic and Algebraic Computation, pages 205-211, New York, 1997. ACM.

[P] A. P. Prudnikov, Yu. A. Brychkov and O. I. Marichev (1990).

Integrals and Series: More Special Functions, Vol. 3,. Gordon and Breach Science Publisher

exception sympy.integrals.meijerint._CoeffExpValueError[source]#

Exception raised by _get_coeff_exp, for internal use only.

sympy.integrals.meijerint._check_antecedents(g1, g2, x)[source]#

Return a condition under which the integral theorem applies.

sympy.integrals.meijerint._check_antecedents_1(g, x, helper=False)[source]#

Return a condition under which the mellin transform of g exists. Any power of x has already been absorbed into the G function, so this is just $$\int_0^\infty g\, dx$$.

See [L, section 5.6.1]. (Note that s=1.)

If helper is True, only check if the MT exists at infinity, i.e. if $$\int_1^\infty g\, dx$$ exists.

sympy.integrals.meijerint._check_antecedents_inversion(g, x)[source]#

Check antecedents for the laplace inversion integral.

sympy.integrals.meijerint._condsimp(cond, first=True)[source]#

Do naive simplifications on cond.

Explanation

Note that this routine is completely ad-hoc, simplification rules being added as need arises rather than following any logical pattern.

Examples

>>> from sympy.integrals.meijerint import _condsimp as simp
>>> from sympy import Or, Eq
>>> from sympy.abc import x, y
>>> simp(Or(x < y, Eq(x, y)))
x <= y

sympy.integrals.meijerint._create_lookup_table(table)[source]#

Add formulae for the function -> meijerg lookup table.

sympy.integrals.meijerint._dummy(name, token, expr, **kwargs)[source]#

Return a dummy. This will return the same dummy if the same token+name is requested more than once, and it is not already in expr. This is for being cache-friendly.

sympy.integrals.meijerint._dummy_(name, token, **kwargs)[source]#

Return a dummy associated to name and token. Same effect as declaring it globally.

sympy.integrals.meijerint._eval_cond(cond)[source]#

Re-evaluate the conditions.

sympy.integrals.meijerint._exponents(expr, x)[source]#

Find the exponents of x (not including zero) in expr.

Examples

>>> from sympy.integrals.meijerint import _exponents
>>> from sympy.abc import x, y
>>> from sympy import sin
>>> _exponents(x, x)
{1}
>>> _exponents(x**2, x)
{2}
>>> _exponents(x**2 + x, x)
{1, 2}
>>> _exponents(x**3*sin(x + x**y) + 1/x, x)
{-1, 1, 3, y}

sympy.integrals.meijerint._find_splitting_points(expr, x)[source]#

Find numbers a such that a linear substitution x -> x + a would (hopefully) simplify expr.

Examples

>>> from sympy.integrals.meijerint import _find_splitting_points as fsp
>>> from sympy import sin
>>> from sympy.abc import x
>>> fsp(x, x)
{0}
>>> fsp((x-1)**3, x)
{1}
>>> fsp(sin(x+3)*x, x)
{-3, 0}

sympy.integrals.meijerint._flip_g(g)[source]#

Turn the G function into one of inverse argument (i.e. G(1/x) -> G’(x))

sympy.integrals.meijerint._functions(expr, x)[source]#

Find the types of functions in expr, to estimate the complexity.

sympy.integrals.meijerint._get_coeff_exp(expr, x)[source]#

When expr is known to be of the form c*x**b, with c and/or b possibly 1, return c, b.

Examples

>>> from sympy.abc import x, a, b
>>> from sympy.integrals.meijerint import _get_coeff_exp
>>> _get_coeff_exp(a*x**b, x)
(a, b)
>>> _get_coeff_exp(x, x)
(1, 1)
>>> _get_coeff_exp(2*x, x)
(2, 1)
>>> _get_coeff_exp(x**3, x)
(1, 3)

sympy.integrals.meijerint._guess_expansion(f, x)[source]#

Try to guess sensible rewritings for integrand f(x).

sympy.integrals.meijerint._inflate_fox_h(g, a)[source]#

Let d denote the integrand in the definition of the G function g. Consider the function H which is defined in the same way, but with integrand d/Gamma(a*s) (contour conventions as usual).

If a is rational, the function H can be written as C*G, for a constant C and a G-function G.

This function returns C, G.

sympy.integrals.meijerint._inflate_g(g, n)[source]#

Return C, h such that h is a G function of argument z**n and g = C*h.

sympy.integrals.meijerint._int0oo(g1, g2, x)[source]#

Express integral from zero to infinity g1*g2 using a G function, assuming the necessary conditions are fulfilled.

Examples

>>> from sympy.integrals.meijerint import _int0oo
>>> from sympy.abc import s, t, m
>>> from sympy import meijerg, S
>>> g1 = meijerg([], [], [-S(1)/2, 0], [], s**2*t/4)
>>> g2 = meijerg([], [], [m/2], [-m/2], t/4)
>>> _int0oo(g1, g2, t)
4*meijerg(((1/2, 0), ()), ((m/2,), (-m/2,)), s**(-2))/s**2

sympy.integrals.meijerint._int0oo_1(g, x)[source]#

Evaluate $$\int_0^\infty g\, dx$$ using G functions, assuming the necessary conditions are fulfilled.

Examples

>>> from sympy.abc import a, b, c, d, x, y
>>> from sympy import meijerg
>>> from sympy.integrals.meijerint import _int0oo_1
>>> _int0oo_1(meijerg([a], [b], [c], [d], x*y), x)
gamma(-a)*gamma(c + 1)/(y*gamma(-d)*gamma(b + 1))

sympy.integrals.meijerint._int_inversion(g, x, t)[source]#

Compute the laplace inversion integral, assuming the formula applies.

sympy.integrals.meijerint._is_analytic(f, x)[source]#

Check if f(x), when expressed using G functions on the positive reals, will in fact agree with the G functions almost everywhere

sympy.integrals.meijerint._meijerint_definite_2(f, x)[source]#

Try to integrate f dx from zero to infinity.

The body of this function computes various ‘simplifications’ f1, f2, … of f (e.g. by calling expand_mul(), trigexpand() - see _guess_expansion) and calls _meijerint_definite_3 with each of these in succession. If _meijerint_definite_3 succeeds with any of the simplified functions, returns this result.

sympy.integrals.meijerint._meijerint_definite_3(f, x)[source]#

Try to integrate f dx from zero to infinity.

This function calls _meijerint_definite_4 to try to compute the integral. If this fails, it tries using linearity.

sympy.integrals.meijerint._meijerint_definite_4(f, x, only_double=False)[source]#

Try to integrate f dx from zero to infinity.

Explanation

This function tries to apply the integration theorems found in literature, i.e. it tries to rewrite f as either one or a product of two G-functions.

The parameter only_double is used internally in the recursive algorithm to disable trying to rewrite f as a single G-function.

sympy.integrals.meijerint._meijerint_indefinite_1(f, x)[source]#

Helper that does not attempt any substitution.

sympy.integrals.meijerint._mul_args(f)[source]#

Return a list L such that Mul(*L) == f.

If f is not a Mul or Pow, L=[f]. If f=g**n for an integer n, L=[g]*n. If f is a Mul, L comes from applying _mul_args to all factors of f.

sympy.integrals.meijerint._mul_as_two_parts(f)[source]#

Find all the ways to split f into a product of two terms. Return None on failure.

Explanation

Although the order is canonical from multiset_partitions, this is not necessarily the best order to process the terms. For example, if the case of len(gs) == 2 is removed and multiset is allowed to sort the terms, some tests fail.

Examples

>>> from sympy.integrals.meijerint import _mul_as_two_parts
>>> from sympy import sin, exp, ordered
>>> from sympy.abc import x
>>> list(ordered(_mul_as_two_parts(x*sin(x)*exp(x))))
[(x, exp(x)*sin(x)), (x*exp(x), sin(x)), (x*sin(x), exp(x))]

sympy.integrals.meijerint._my_principal_branch(expr, period, full_pb=False)[source]#

Bring expr nearer to its principal branch by removing superfluous factors. This function does not guarantee to yield the principal branch, to avoid introducing opaque principal_branch() objects, unless full_pb=True.

sympy.integrals.meijerint._mytype(f, x)[source]#

Create a hashable entity describing the type of f.

sympy.integrals.meijerint._rewrite1(f, x, recursive=True)[source]#

Try to rewrite f using a (sum of) single G functions with argument a*x**b. Return fac, po, g such that f = fac*po*g, fac is independent of x. and po = x**s. Here g is a result from _rewrite_single. Return None on failure.

sympy.integrals.meijerint._rewrite2(f, x)[source]#

Try to rewrite f as a product of two G functions of arguments a*x**b. Return fac, po, g1, g2 such that f = fac*po*g1*g2, where fac is independent of x and po is x**s. Here g1 and g2 are results of _rewrite_single. Returns None on failure.

sympy.integrals.meijerint._rewrite_inversion(fac, po, g, x)[source]#

Absorb po == x**s into g.

sympy.integrals.meijerint._rewrite_saxena(fac, po, g1, g2, x, full_pb=False)[source]#

Rewrite the integral fac*po*g1*g2 from 0 to oo in terms of G functions with argument c*x.

Explanation

Return C, f1, f2 such that integral C f1 f2 from 0 to infinity equals integral fac po, g1, g2 from 0 to infinity.

Examples

>>> from sympy.integrals.meijerint import _rewrite_saxena
>>> from sympy.abc import s, t, m
>>> from sympy import meijerg
>>> g1 = meijerg([], [], [0], [], s*t)
>>> g2 = meijerg([], [], [m/2], [-m/2], t**2/4)
>>> r = _rewrite_saxena(1, t**0, g1, g2, t)
>>> r[0]
s/(4*sqrt(pi))
>>> r[1]
meijerg(((), ()), ((-1/2, 0), ()), s**2*t/4)
>>> r[2]
meijerg(((), ()), ((m/2,), (-m/2,)), t/4)

sympy.integrals.meijerint._rewrite_saxena_1(fac, po, g, x)[source]#

Rewrite the integral fac*po*g dx, from zero to infinity, as integral fac*G, where G has argument a*x. Note po=x**s. Return fac, G.

sympy.integrals.meijerint._rewrite_single(f, x, recursive=True)[source]#

Try to rewrite f as a sum of single G functions of the form C*x**s*G(a*x**b), where b is a rational number and C is independent of x. We guarantee that result.argument.as_coeff_mul(x) returns (a, (x**b,)) or (a, ()). Returns a list of tuples (C, s, G) and a condition cond. Returns None on failure.

sympy.integrals.meijerint._split_mul(f, x)[source]#

Split expression f into fac, po, g, where fac is a constant factor, po = x**s for some s independent of s, and g is “the rest”.

Examples

>>> from sympy.integrals.meijerint import _split_mul
>>> from sympy import sin
>>> from sympy.abc import s, x
>>> _split_mul((3*x)**s*sin(x**2)*x, x)
(3**s, x*x**s, sin(x**2))

sympy.integrals.meijerint.meijerint_definite(f, x, a, b)[source]#

Integrate f over the interval [a, b], by rewriting it as a product of two G functions, or as a single G function.

Return res, cond, where cond are convergence conditions.

Examples

>>> from sympy.integrals.meijerint import meijerint_definite
>>> from sympy import exp, oo
>>> from sympy.abc import x
>>> meijerint_definite(exp(-x**2), x, -oo, oo)
(sqrt(pi), True)

This function is implemented as a succession of functions meijerint_definite, _meijerint_definite_2, _meijerint_definite_3, _meijerint_definite_4. Each function in the list calls the next one (presumably) several times. This means that calling meijerint_definite can be very costly.

sympy.integrals.meijerint.meijerint_indefinite(f, x)[source]#

Compute an indefinite integral of f by rewriting it as a G function.

Examples

>>> from sympy.integrals.meijerint import meijerint_indefinite
>>> from sympy import sin
>>> from sympy.abc import x
>>> meijerint_indefinite(sin(x), x)
-cos(x)

sympy.integrals.meijerint.meijerint_inversion(f, x, t)[source]#

Compute the inverse laplace transform $$\int_{c+i\infty}^{c-i\infty} f(x) e^{tx}\, dx$$, for real c larger than the real part of all singularities of f.

Note that t is always assumed real and positive.

Return None if the integral does not exist or could not be evaluated.

Examples

>>> from sympy.abc import x, t
>>> from sympy.integrals.meijerint import meijerint_inversion
>>> meijerint_inversion(1/x, x, t)
Heaviside(t)


text:
Context.—

Postmortem evaluation for neurodegenerative disease is expensive in time and materials. These challenges can be met by implementing simpler sampling protocols while preserving anatomic relations.

Objective.—

To determine the diagnostic effectiveness and cost-effectiveness of a simplified brain blocking protocol compared with the standard blocking protocol used in our Alzheimer's Disease Research Center (ADRC).

Design.—

We prospectively compared the neuropathologic diagnoses established from our standard 19-cassette/19 brain sites ADRC protocol to a simplified 6-cassette/12 brain sites protocol in 52 consecutive cases. The simplified protocol generated 14 slides for comparison to 52 slides from our standard protocol.

Results.—

Compared with the ADRC protocol the simplified protocol produced Alzheimer Disease Neuropathologic Changes probability scores that were the same in 50 of 52 cases (r = 0.99). Staging for Lewy pathology was equivalent in 45 of 52 (r = 0.98), scoring for cerebral amyloid angiopathy was equivalent in 48 of 52 (r = 0.97), and grading for arteriolosclerosis was the same in 45 of 52 cases (r = 0.92). Progressive supranuclear palsy (n = 4), multiple system atrophy (n = 2), and corticobasal degeneration (n = 1) could be diagnosed by either protocol independently. The estimated savings per case was 72% or $1744.89 ($2436.37 [ADRC] versus $691.48 [simplified]). Conclusions.— The diagnosis of neurodegenerative disease at autopsy can be done accurately with a less expensive, simplified protocol. Our protocol is similar to those of previously published approaches, but it has a simpler organization scheme. This method should be valuable to institutions where autopsy cost considerations may be important. Limited funding for autopsies generates significant monetary concerns for any center with a busy autopsy service. In the period 2017–2018 at the Massachusetts General Hospital, Boston Massachusetts, we performed 852 autopsies, 772 of which included examination of the brain. Although 205 of these cases were part of our Massachusetts Alzheimer's Disease Research Center (MADRC), Boston Massachusetts, and received associated funding for processing of tissue, the other 567 received no outside support. Because Alzheimer disease (AD) is highly prevalent in the United States (estimated 5.7 million cases in 20181 ) and its burden increases significantly with age (3% of those aged 65–74 years, 17% of those aged 75–84 years, and 32% of those 85 years and older1 ), the ability to routinely screen for AD and other neuropathologic diseases at the time of autopsy is important from a diagnostic and public health standpoint. In a climate of increasing health costs, combined with a lack of reimbursement for autopsies, economic assessment of diagnostic methods is essential. In 2018 we attempted to implement a condensed brain blocking protocol that results in considerable cost savings.2 We found this protocol to be cumbersome to prepare because it required sampling of 20 anatomic brain sites, with each piece assigned to a specifically designated place in 1 of 5 cassettes (4 brain sites per cassette). The small size of each piece limited examination of the anatomy. For example, the anatomic relationship between the hippocampus, parahippocampal gyrus, and occipital-temporal gyrus was not preserved, an important factor when performing Braak staging for neurofibrillary tangles. We also noted redundancies in the sampling (ie, bilateral thalamic sampling), which may be unnecessary for the diagnosis of neurodegenerative disease. To avoid these issues while preserving the ability to reduce cost, we developed a simplified blocking protocol that required sampling 12 brain sites, with only 2 pieces of tissue placed in each cassette for a total of 6 cassettes. We found this method much easier to prepare, the blocking more intuitive to remember, and the anatomic relationships within each tissue piece better preserved. We validated the performance of this protocol by comparing it to our standard 19-cassette/19 brain sites MADRC protocol in 52 consecutive cases submitted to MADRC. We studied 52 consecutive submissions to our MADRC during the period of approximately 8 months with antemortem diagnoses of suspected AD (n = 18), parkinsonian-related disorder (n = 10), primary progressive aphasia (n = 6), frontotemporal dementia (n = 4), mixed dementia (n = 4), corticobasal syndrome (n = 3), cerebral amyloid angiopathy (n = 3), mild cognitive impairment (n = 3), multiple system atrophy (n = 2), normal pressure hydrocephalus (n = 1), and cases submitted as controls (n = 5; Table 3). We excluded from this study cases submitted for evaluation of amyotrophic lateral sclerosis, Huntington disease, and spinal muscular atrophy because we use a separate blocking protocol at our Alzheimer's Disease Research Center (ADRC). Consent for autopsy, including consent for the processing of brain for research purposes, was obtained from the legal health care proxy in each case. Each brain was processed according to our standard MADRC protocol, with the brain divided at the midline at the time of autopsy, with one-half cut into coronal slabs for freezing at −80°C and the other half fixed in 10% buffered formalin. Following at least 1 week of formalin fixation, the fixed hemisphere was sectioned in the coronal plane at standard landmarks as previously described.3 Twelve brain sites were sampled from the formalin-fixed side and placed in 6 cassettes (hereafter called simplified protocol; Figure 1 and Table 1). Nineteen additional brain sites were sampled from the formalin-fixed side and placed in 19 cassettes for further processing, staining, and microscopic analysis, according to our standard ADRC protocol (Figure 2 and Table 2). Histologic sections corresponding to blocks in Figure 1 (1 through 6) are seen in Figure 3 (A through F). Figure 1 Cassettes used for the simplified protocol (numbered 1–6). 1, Cerebellum; basal ganglia; 2, hippocampus at lateral geniculate nucleus (LGN); superior frontal cortex (Brodmann Area [BA] 8, 9); 3, midbrain; cingulate gyrus (BA24); 4, medulla oblongata; pons; 5, calcarine cortex (BA17); superior frontal cortex (BA7); 6, thalamus; temporal pole (BA38). See Table 1 for staining protocol. Figure 1 Cassettes used for the simplified protocol (numbered 1–6). 1, Cerebellum; basal ganglia; 2, hippocampus at lateral geniculate nucleus (LGN); superior frontal cortex (Brodmann Area [BA] 8, 9); 3, midbrain; cingulate gyrus (BA24); 4, medulla oblongata; pons; 5, calcarine cortex (BA17); superior frontal cortex (BA7); 6, thalamus; temporal pole (BA38). See Table 1 for staining protocol. Table 1 Simplified Staining Protocola Figure 2 Cassettes used for the MADRC protocol (numbered 1–19). 1, Superior frontal cortex (Brodmann Area [BA] 8, 9); 2, primary motor/sensory strip (BA 3, 4); 3, superior parietal cortex (BA7); 4, calcarine cortex (BA17); 5, hippocampus at lateral geniculate nucleus (LGN); 6, caudate, putamen, nucleus accumbens; 7, globus pallidus, putamen, substantia innominate; 8, amygdala and entorhinal cortex; 9, thalamus (level of centrum medianum); 10, midbrain; 11, upper pons (level of locus ceruleus); 12, lower pons (at inferior border of cranial nerve V); 13, medulla oblongata; 14, cerebellum (with dentate nucleus); 15, temporal pole (BA38); 16, cingulate gyrus (BA24); 17, thalamus with anterior nucleus, ventral anterior, ventral lateral, and subthalamic nucleus; 18, inferior frontal cortex (BA 10, 11, 12); 19, inferior parietal cortex (BA 39, 40). See Table 2 for staining protocol. Figure 2 Cassettes used for the MADRC protocol (numbered 1–19). 1, Superior frontal cortex (Brodmann Area [BA] 8, 9); 2, primary motor/sensory strip (BA 3, 4); 3, superior parietal cortex (BA7); 4, calcarine cortex (BA17); 5, hippocampus at lateral geniculate nucleus (LGN); 6, caudate, putamen, nucleus accumbens; 7, globus pallidus, putamen, substantia innominate; 8, amygdala and entorhinal cortex; 9, thalamus (level of centrum medianum); 10, midbrain; 11, upper pons (level of locus ceruleus); 12, lower pons (at inferior border of cranial nerve V); 13, medulla oblongata; 14, cerebellum (with dentate nucleus); 15, temporal pole (BA38); 16, cingulate gyrus (BA24); 17, thalamus with anterior nucleus, ventral anterior, ventral lateral, and subthalamic nucleus; 18, inferior frontal cortex (BA 10, 11, 12); 19, inferior parietal cortex (BA 39, 40). See Table 2 for staining protocol. Table 2 Massachusetts Alzheimer's Disease Research Center Staining Protocola Table 3 Case Details Figure 3 Microscopy from simplified protocol. A, Cerebellum; basal ganglia. B, Hippocampus at lateral geniculate nucleus (LGN); superior frontal cortex (Brodmann Area [BA] 8, 9). C, Midbrain; cingulate gyrus (BA24). D, Medulla oblongata; pons. E, Calcarine cortex (BA17); superior parietal cortex (BA7). F, Thalamus; temporal pole (BA38) (Luxol Fast Blue and hematoxylin-eosin, original magnification ×1). Figure 3 Microscopy from simplified protocol. A, Cerebellum; basal ganglia. B, Hippocampus at lateral geniculate nucleus (LGN); superior frontal cortex (Brodmann Area [BA] 8, 9). C, Midbrain; cingulate gyrus (BA24). D, Medulla oblongata; pons. E, Calcarine cortex (BA17); superior parietal cortex (BA7). F, Thalamus; temporal pole (BA38) (Luxol Fast Blue and hematoxylin-eosin, original magnification ×1). Tissues for the simplified protocol were processed in the clinical laboratory at Massachusetts General Hospital on a Leica Peloris tissue processor (Leica Biosystems Inc, Buffalo Grove, Illinois) and were embedded in paraffin. Sections of 5 μm sections were stained with Luxol Fast Blue and hematoxylin-eosin (L/H&E) for routine assessment. Eight μm sections from selected blocks were immunostained for hyperphosphorylated tau (Polyclonal Rabbit Anti-Human Tau, Dako, Glostrup, Denmark) at a titration of 1:150 following antigen retrieval with Enzyme 1, β-amyloid (rabbit monoclonal, clone RBT-A4, BioSB, Santa Barbara, California) ready-to-use following HIER for 2.5 minutes with ER2, and α-synuclein (mouse monoclonal, clone BSB-114, BioSB) ready-to-use following HIER for 20 minutes with ER1 using a Leica Bond III automated stainer (Leica Biosystems, Wetzlar, Germany) according to the manufacturer's instructions. The simplified protocol resulted in a total of 6 L/H&E, 2 tau, 3 β-amyloid, and 3 α-synuclein stained slides (14 slides in total; Table 1). Tissues for the MADRC protocol were handled in our research lab. The formalin-fixed blocks were treated with 88% formic acid for 60 minutes followed by processing on a Thermo Scientific Excelsior ES tissue processor (Thermo Fisher Scientific, Waltham, Massachusetts) and embedded in paraffin. Sections 7 μm thick were stained with L/H&E for routine assessment. Staining for Bielschowsky silver stain was performed on sections from selected blocks (Table 2). Immunohistochemistry for hyperphosphorylated tau (Polyclonal Rabbit Anti-Human Tau) at a titration of 1:6000, β-amyloid (Monoclonal mouse Anti-Human Beta-Amyloid clone 6F/3D, Dako) at a titration of 1:600, α-synuclein (Mouse anti-α-Synuclein LB509, Life Technologies Corp, Frederick, Maryland) at a titration of 1:200, and p-TDP-43 (Cosmo Bio, Tokyo, Japan) at a titration of 1:3000 were performed on 7-μm–thick sections from select blocks and processed on a Leica Bond RX automated stainer (Leica Biosystems, Wetzlar, Germany) according to the manufacturer's instructions. The MADRC protocol resulted in 19 L/H&E, 8 Bielschowsky silver, 7 tau, 6 α-synuclein, 10 β-amyloid, and 2 p-TDP-43 stained slides (52 slides in total; Table 2). All cases were prospectively analyzed by 1 neuropathology fellow and a neuropathology attending pathologist. The MADRC and simplified protocols were evaluated concurrently but scored independently for evidence of neuropathologic disease. The extent of Alzheimer disease neuropathologic changes (ADNCs) was established using the 2012 NIAA guidelines.4 Staining for β-amyloid was used for Thal phase (scale, 0–5).5 Bielschowsky silver stain and hyperphosphorylated tau stains were used for Braak and Braak neurofibrillary tangles stage (scale, none to VI) and CERAD neuritic plaque score (scale, none to frequent) in the ADRC protocol.68 Only the hyperphosphorylated tau stain was used for neurofibrillary tangles staging and CERAD neuritic plaque score in the simplified protocol. The Thal, neurofibrillary tangles stage, and CERAD scores were used to generate an “ABC” score (scale, 0–3), which in turn was used to determine the probability of ADNC (not, low, intermediate, and high). The extent of Lewy body disease was similarly determined and scored according to Braak staging (scale, 0–6) for each protocol.911 Comparisons for cerebrovascular disease were also performed. The extent of cerebral amyloid angiopathy was scored on a 5-point scale (scale, 0–4; 0 = none, 1 = mild, 2 = moderate, 3 = severe, 4 = very severe) as described previously.12,13 Arteriolosclerosis was scored as none, mild, moderate, or severe. Microinfarcts were scored as present or not. Comparisons for macroinfarcts were not done because these require targeted sampling at the time of grossing which deviates from either protocol. Evaluation for presence of abnormal TDP-43 staining was done in the MADRC protocol. From the cases that were positive for TDP-43 with the MADRC protocol, sections from the simplified protocol-hippocampal block were stained for TDP-43 at the research lab, for a retrospective comparison. The cases included 28 men (28 of 52; 53.8%) and 24 women (24 of 52; 46.2%) with a mean age of 78 years (range, 52–102 years), mean pre–formalin fixation brain weight of 1170 g (range, 820–1520 g), and mean postmortem interval time of 22 hours (range, 4–138 hours; Table 3). The anatomic relationships were well preserved with the simplified protocol (Figure 1). The cassette with the hippocampus and superior frontal gyrus allowed for preservation of the relationship of Ammon horn with the parahippocampal and occipital-temporal gyri as well as in some instances (depending on the degree of atrophy) a portion of the inferior temporal gyrus (Figure 1, cassette 2). The cassette with the cerebellum and basal ganglia included the cerebellar cortex, deep cerebellar white matter, and a portion of the dentate nucleus (Figure 1, cassette 1). This cassette also allowed for representation of the caudate, internal capsule, and putamen, and in some instances the globus pallidus. The cassette with the temporal cortex included thalamus and subthalamic nucleus (Figure 1, cassette 6). This area was helpful in the evaluation of progressive supranuclear palsy and corticobasal degeneration. Taking a smaller section of cingulate gyrus allowed for the entire axial half of the midbrain to be sampled in the same cassette (Figure 1, cassette 3). The cassette with the calcarine sulcus included the primary and association visual cortices (allowing for a determination of Braak V versus Braak VI) as well as the separately sampled superior parietal cortex (Figure 1, cassette 5). The overall probability scores for ADNC were the same for the ADRC and simplified protocols in 50 of 52 cases (98% accuracy, r = 0.99), with equivalent scores of “not” (n = 9), “low” (n = 13), “intermediate” (n = 3), and “high” (n = 25; Table 4). The 2 cases that differed were both categorized as “low” with the simplified protocol and “intermediate” with the ADRC protocol. There were no differences in the “not” and “high” categories across protocols (Table 4). Examples of Thal 5, Braak VI, and CERAD frequent via immunohistochemistry from the simplified protocol are shown in Figure 4, A through C. Table 4 Neuropathology Results, Massachusetts Alzheimer's Disease Research Center (MADRC; 19) Versus Simplified (6)a Figure 4 Immunohistochemistry examples from the simplified protocol. A, β-Amyloid with diffuse plaques involving cerebellum, example of Thal 5, “A3.” B, p-Tau with neurofibrillary tangles involving primary visual cortex, example of Braak VI, “B3.” C, p-Tau showing frequent neuritic plaques (arrows point to a few of these for emphasis) involving neocortex, CERAD frequent, “C3.” D, α-Synuclein with neocortical Lewy bodies in a case of Lewy body disease. E, α-Synuclein with a single Lewy body involving Cornu Ammonis subfield 1 (CA1) of the hippocampus, a location associated with Alzheimer disease. F, α-Synuclein showing cerebellar glial cytoplasmic inclusions in a case of multiple system atrophy. G, p-Tau showing a tufted astrocyte in progressive supranuclear palsy. H, p-Tau showing globose tangle in a case of corticobasal degeneration. I, β-Amyloid with cortical and leptomeningeal involvement by cerebral amyloid angiopathy (CAA). J, β-Amyloid with dyshoric CAA. K, p-TDP-43 showing neuronal cytoplasmic inlcusions (NCIs) in a case of frontotemporal lobar degeneration with TDP-43. L, p-TDP-43 demonstrating NCIs in the dentate gyrus of the hippocampus, in a location associated with Alzheimer disease (original magnifications ×200 [A, C, and F], ×400 [B, D, E, G, H, J, K, L], and ×40 [I]). Figure 4 Immunohistochemistry examples from the simplified protocol. A, β-Amyloid with diffuse plaques involving cerebellum, example of Thal 5, “A3.” B, p-Tau with neurofibrillary tangles involving primary visual cortex, example of Braak VI, “B3.” C, p-Tau showing frequent neuritic plaques (arrows point to a few of these for emphasis) involving neocortex, CERAD frequent, “C3.” D, α-Synuclein with neocortical Lewy bodies in a case of Lewy body disease. E, α-Synuclein with a single Lewy body involving Cornu Ammonis subfield 1 (CA1) of the hippocampus, a location associated with Alzheimer disease. F, α-Synuclein showing cerebellar glial cytoplasmic inclusions in a case of multiple system atrophy. G, p-Tau showing a tufted astrocyte in progressive supranuclear palsy. H, p-Tau showing globose tangle in a case of corticobasal degeneration. I, β-Amyloid with cortical and leptomeningeal involvement by cerebral amyloid angiopathy (CAA). J, β-Amyloid with dyshoric CAA. K, p-TDP-43 showing neuronal cytoplasmic inlcusions (NCIs) in a case of frontotemporal lobar degeneration with TDP-43. L, p-TDP-43 demonstrating NCIs in the dentate gyrus of the hippocampus, in a location associated with Alzheimer disease (original magnifications ×200 [A, C, and F], ×400 [B, D, E, G, H, J, K, L], and ×40 [I]). Assessment of Lewy pathology was the same for the MADRC and simplified protocols in 45 of 52 cases (87% accuracy), with equivalent scores of “none” (n = 30), Braak stages 3 (n = 1), 4 (n = 1), and 5 (n = 9), as shown in Figure 4, D; and 6 (n = 1), and “AD-associated” (n = 3), as shown in Figure 4, E. The 6 cases that differed were Braak and Braak differences of 5 versus 6 (n = 1), 4 versus 5 (n = 3), 3 versus 5 (n = 1), and “none” versus “AD associated” (n = 2). In each of the cases with differences in Braak or Lewy body scores, the simplified protocol had a lower score, but for cases assigned a Braak stage, the correlation was high (r = 0.98; Table 4). Multiple system atrophy (n = 2; Figure 4, F), progressive supranuclear palsy (n = 4; Figure 4, G), and corticobasal degeneration (n = 1; Figure 4, H) were diagnosed with either protocol. Scoring for cerebral amyloid angiopathy (Figure 4, I and J) was equivalent in 48 of 52 cases (92% accuracy, r = 0.97) with scores of “none” (n = 28), grade 1 (n = 2), grade 2 (n = 11), grade 3 (n = 6), and grade 4 (n = 1). Four cases differed by 1 degree only, grade 1 versus grade 2 (n = 3) and grade 2 versus grade 3 (n = 1; Table 4). Grading for arteriolosclerosis was the same in 45 of 52 cases (r = 0.92) with equivalent scores of “none” (n = 4), mild (n = 16), moderate (n = 16), and severe (n = 3). Seven cases differed by 1 degree each (Table 4), with the MADRC protocol showing a higher degree. Microinfarcts were present in 5 of 52 cases with the MADRC protocol, and in 3 of 52 cases with the simplified protocol (Table 4). TDP-43 positivity was present in 15 of 52 cases (29%) stained in the MADRC protocol. Three cases demonstrated neuropathologic changes of FTLD-TDP-43, as shown in Figure 4, K. The other 12 had inclusions and neurites in the entorhinal cortex, amygdala, and/or hippocampus in a pattern associated with ADNC (Figure 4, L).14 Staining for TDP-43 was performed on the hippocampus/superior frontal block from the simplified protocol for these 15 cases. Of these, 13 (87%) had TDP-43 inclusions and neurites. The 2 cases that were negative in the simplified protocol either had only focal staining for TDP-43 in the amygdala and anterior entorhinal cortex (an area not represented with the simplified protocol) or had only focal TDP-43 positivity in pyramidal neurons of CA1 with the MADRC protocol. The estimated cost for carrying out the MADRC protocol in our clinical laboratory was$2436.37, whereas that for the simplified protocol was only $691.48 (Table 5). The$1744.89 difference in cost represents a 72% reduction in expenses. Most of the savings were in the reduced costs for immunohistochemistry with the simplified protocol. The number of slides to review was decreased by 73% (14 versus 52 slides), saving significant diagnostic time compared with the MADRC protocol, which is both a diagnostic and research protocol.

Table 5

Analysis of Estimated Cost Savingsa

We have developed a cost-effective, simplified protocol, similar to the condensed technique described by Flanagan et al2  but modified with a simplified organizational scheme and larger individual pieces with better preservation of the neuroanatomy. We have validated our technique in a prospective cohort of cases submitted to the MADRC by comparing it with the combined diagnostic and research protocol as the “gold standard.” Our technique also has excellent correlation with the MADRC protocol for cerebrovascular disease, especially for cerebral amyloid angiopathy and arteriolosclerosis. This is most likely due to the larger size samples compared with previously published protocols because each of our cortical blocks included leptomeninges and subcortical white matter. Microinfarcts were detected at a slightly reduced rate compared with the MADRC protocol, likely attributable to the smaller overall volume of sampling, a limitation also present with the Flanagan method.

This protocol is optimized for the diagnosis of neurodegenerative disease. It is 98% accurate in the assessment of ADNC in our series and 100% accurate when the classification is “high” probability of AD or effectively “not” AD. The accuracy of an “intermediate” or “low” probability classification is 89%. This reduction in accuracy is of less diagnostic relevance given that an “intermediate” or “low” ADNC score will not definitively rule in or rule out the diagnosis of AD.

The simplified protocol is also very effective in detecting the presence of Lewy pathology, with an overall accuracy of 87%, with some decreased sensitivity attributed to reduced overall sampling compared with the MADRC protocol. All cases that were submitted with a clinical suspicion of Lewy pathology, as well as those that had evidence of Lewy body disease with the MADRC protocol, also had evidence of Lewy body disease with the simplified protocol. Hence, no false negatives were detected.

Although we did not prospectively test for TDP-43 with the simplified protocol, the hippocampal/superior frontal block should be stained for TDP-43 to screen for hippocampal sclerosis of aging, FTLD-TDP-43, and the TDP-43 staining associated with ADNC. In our retrospective review, 87% of the cases that showed positivity for TDP-43 in the MADRC protocol also showed positivity with the simplified protocol. Those that did not stain with the simplified protocol were cases with only focal positive staining with the MADRC protocol. In clinical practice, this stain should be performed on block 2 of the simplified protocol when the other immunohistochemical stains (β-amyloid, hyperphosphorylated tau, α-synuclein) are used.

Less common neurodegenerative disorders were also effectively diagnosed with the simplified protocol. Our cohort of cases included a number of cases with progressive supranuclear palsy and multiple system atrophy, and 1 case of corticobasal degeneration. The simplified blocking scheme allowed sufficient tissue for additional immunohistochemistry to be done when needed, with minimal additional cost. For example, for the case of multiple system atrophy we added α-synuclein to block 1 of the simplified protocol to identify glial inclusions in the cerebellar white matter, internal capsule, and putamen. Similarly, if there is concern for cerebral amyloid angiopathy, the protocol could be modified to add a β-amyloid stain to block 5 (superior parietal, occipital cortex) to increase diagnostic yield because these regions are often involved.

The potential savings when using the simplified protocol are substantial. The 72% potential reduction in cost is significant when running a busy autopsy service that is generally not directly reimbursed, or in cost-restrictive settings, such as small rural academic centers, or a medical examiner's office. This reduction in cost is similar to a previously published alternative way of reducing the amount of tissue processed for microscopy.2  The considerable time saved by the pathologist looking at the case, having reduced the number of slides needed to make the diagnosis, is also substantial, as is the reduced technician time required in embedding and cutting sections.

Several of the authors have implemented this simplified protocol in clinical autopsy practice for older individuals where a clinical history of dementia is lacking or unclear but may potentially be of concern. The L/H&E stains are reviewed from the 6 blocks for plaques, tangles, Lewy bodies, and other pathology, and then a decision as to whether to order the immunohistochemistry is made. In many cases the immunohistochemistry is not needed and additional savings in time and money are gained.

A limitation of our study is that although we examined all cases prospectively, in most instances we were not blinded when comparing the scores across protocols. However, all slides were reviewed by a neuropathology attending with considerable experience in the field of neurodegenerative disease. Another limitation was that we did not sample the amygdala and anterior transentorhinal region in the simplified protocol. This blocking scheme would miss p-tau positive neurofibrillary tangles and/or p-TDP-43 positive inclusions isolated to those regions, as well as amygdala-only Lewy body disease. However, these scenarios may be of less clinical interest when these regions are involved in isolation. A modification of the protocol could potentially substitute the amygdala for the thalamus because the thalamus was often helpful in only a minority of cases.

A third limitation was that the stains were done in separate facilities. This, however, gave us a chance to carefully look at the differences and see where improvements in staining protocols could be made.

Although we excluded cases of Huntington disease, amyotrophic lateral sclerosis, and spinal muscular atrophy from our study, because we use different blocking protocols, a modification of the simplified protocol could be used to examine these cases. For example, in amyotrophic lateral sclerosis, a seventh block could be added to the 6-block protocol. This cassette would contain a section from the precentral gyrus and a single axial section each of cervical, thoracic, and lumbar cord. An L/H&E, glial fibrillary acid protein, ubiquitin, and TDP-43 could be done on this block, providing evidence of upper and lower motor neuron disease as well as the inclusions seen in this disease.

In summary, in any practice setting with a large volume of autopsies in older patients, the need for an economical way to rapidly assess for neurodegenerative disease is essential. This is especially important in cases where the clinical history of dementia is uncertain or lacking. This simplified protocol is sufficient to establish the diagnosis of ADNC, Lewy body pathology, cerebral amyloid angiopathy and arteriolosclerosis, as well as more rare neurodegenerative pathology. In addition, the sample sizes are volumetrically adequate to evaluate anatomic relationships necessary for the evaluation of these diseases.

1.
Alzheimer's Association.
2018 Alzheimer's disease facts and figures
.
Alzheimers Dement
.
2018
;
14
(3)
:
367
429
.
2.
Flanagan
ME,
Marshall
DA,
Shofer
JB,
et al
Performance of a condensed protocol that reduces effort and cost of NIA-AA guidelines for neuropathologic assessment of Alzheimer disease
.
J Neuropathol Exp Neurol
.
2017
;
76
(1)
:
39
43
.
3.
Vonsattel
JP,
Aizawa
H,
Ge
P,
et al
An improved approach to prepare human brains for research
.
J Neuropathol Exp Neurol
.
1995
;
54
(1)
:
42
56
.
4.
Hyman
BT,
Phelps
CH,
Beach
TG,
et al
National Institute on Aging-Alzheimer's Association guidelines for the neuropathologic assessment of Alzheimer's disease
.
Alzheimers Dement
.
2012
;
8
(1)
:
1
13
.
5.
Thal
DR,
Rüb
U,
Orantes
M,
et al
Phases of A beta-deposition in the human brain and its relevance for the development of AD
.
Neurology
.
2002
;
58
(12)
:
1791
1800
.
6.
Braak
H,
Braak
E.
Neuropathological staging of Alzheimer-related changes
.
Acta Neuropathol
.
1991
;
82
(4)
:
239
259
.
7.
Mirra
SS,
Heyman
A,
McKeel
D,
et al
The Consortium to Establish a Registry for Alzheimer Disease (CERAD), part II: standardization of the neuropathologic assessment of Alzheimer's disease
.
Neurology
.
1991
;
41
(4)
:
479
486
.
8.
Montine
TJ,
Phelps
CH,
Beach
TG,
et al
National Institute on Aging-Alzheimer's Association guidelines for the neuropathologic assessment of Alzheimer's disease: a practical approach
.
Acta Neuropathol
.
2012
;
123
(1)
:
1
11
.
9.
Braak
H,
DelTredici
,
Rub
U,
et al
Staging of brain pathology related to sporadic Parkinson's disease
.
Neurobiol Aging
.
2003
;
24
(2)
:
197
211
.
10.
Del Tredici
K,
Rüb
U,
De Vos
RA,
et al
Where does Parkinson disease pathology begin in the brain?
J Neuropathol Exp Neurol
.
2002
;
61
(5)
:
413
426
.
11
McKeith
IG,
Dickson
DW,
Lowe
J,
et al
Diagnosis and management of dementia with Lewy bodies: third report of the DLB Consortium
.
Neurology
.
2005
;
65
(12)
:
1863
1872
.
12.
Greenberg
SM,
Vonsattel
JP.
Diagnosis of cerebral amyloid angiopathy: sensitivity and specificity of cortical biopsy
.
Stroke
.
1997
;
28
(7)
:
1418
1422
.
13.
Arvanitakis
Z,
Leurgans
SE,
Wang
Z.
Cerebral amyloid angiopathy pathology and cognitive domains in older persons
.
Ann Neurol
.
2011
;
69
(2)
:
320
327
.
14.
Josephs
KA,
Melissa
E.
Murray
ME,
et al
Staging TDP-43 pathology in Alzheimer's disease
.
Acta Neuropathol
.
2014
;
127
(3)
:
441
450
.

## Author notes

This work was supported by grant award No. P50 AG005134.

## Competing Interests

The authors have no relevant financial interest in the products or companies described in this article.

Presented in part as a poster presentation at the annual meeting of the American Association for Neuropathologists; June 8, 2019; Atlanta, Georgia.

This work was supported by the Department of Pathology, Brooke Army Medical Center, San Antonio, Texas. The view(s) expressed herein are those of the author(s) and do not reflect the official policy or position of Brooke Army Medical Center, the US Army Medical Department, the US Army Office of the Surgeon General, the Department of the Air Force, the Department of the Army, the Department of Defense, or the US government.


text:
PhysicsBoost! Focus on Astronomy

Physics is sometimes called the philosophy of science.  It's abstract and mathematical and hard for a lot of people. But the concepts behind physics (force, acceleration, friction, etc.) are the forces that make our everyday world work (or fall in on itself).  And, Astronomy is one of the facets of physics that most people are enthralled by.

This course will teach students to use math, physics, and hands-on activities to learn how things work.  They will experiement with simple machines like pulleys, put forces to work in mechanisms like "air cannons," and explore space. As they create and explore, students will build their math skills and their understanding of how and why to use formulae as they apply concepts and ideas with their hands and basic materials.

This class is taking the new Mars rover as a starting point to learn about space: what does the solar system look like?  How do we study it?  How does the sun work?  What does it look like (no, you should never look straight at the sun, but with the right equipment, you can study it safely!)? And, once we get a grip on our solar system, what's beyond?  What are the theories?  And what do current explorations look like?

Sign up using the links below.  You can also sign up for aftercare (3pm-6pm). You can sign up for the week ($100) or for a few days (you can let us know which days on your form, or closer to the date of the camp). Grades: 4th 5th 6th 7th 8th Cost:$350 (Scholarships available)
Time:
9am - 3pm, M-F
Sessions:
Monday, August 9, 2021 to Friday, August 13, 2021