Test which distribution has a "longer tail" I have measured two non-negative random variables, A and B. Their true underlying probabilities are unknown, however, it may be assumed that the probabilities are largest at zero and monotonically decrease for larger values. Most certainly, those values also have an upper bound and I have a guess for it, but that guess is not very good.
I would like to test if the "tail" of $P(A)$ "goes further" than the tail of $P(B)$. It looks like it does, but maybe that is by chance? What metrics could I consider? I have tried to check the mean, but it seems comparable for both variables. 
 A: The essential features of this question are:


*

*It does not make strong distributional assumptions, lending it a non-parametric flavor.

*It concerns only tail behavior, not the entire distribution.
With some diffidence--because I have not studied my proposal theoretically to fully understand its performance--I will outline an approach that might be practicable.  It borrows from the concepts behind the Kolmogorov-Smirnov test, familiar rank-based non-parametric tests, and exploratory data analysis methods.
Let's begin by visualizing the problem.  We may plot the empirical distribution functions of the datasets on common axes to compare them:

The black curve shows dataset $A$ (here with $m=50$ values) and the red curve shows dataset $B$ (here with $n=100$ values).  The height of a curve at a value $x$ shows the proportion of the dataset with values less than or equal to $x.$
This is a situation where data in the upper half of $A$ consistently exceed the data in the upper half of $B.$  We can see that because, scanning from left to right (low values to high values), the curves last cross around a height of $0.5$ and after that, the curve for $A$ (black) remains to the right of -- that is, at higher values than -- the curve for $B$ (red).  That's evidence for a heavier right tail in the distribution from which data $A$ are drawn.
We need a test statistic. It must be a way of somehow quantifying whether and by how much $A$ has a "heavier right tail" than $B.$  My proposal is this:


*

*Combine the two datasets into a dataset of $n+m$ values.

*Rank them: this assigns the value $n+m$ to the highest, $n+m-1$ to the next highest, and so on down to the value $1$ for the lowest.

*Weight the ranks as follows:


*

*Divide the ranks for $A$ by $m$ and the ranks for $B$ by $n.$

*Negate the results for $B.$


*Accumulate these values (in a cumulative sum), beginning with the largest rank and moving on down.

*Optionally, normalize the cumulative sum by multiplying all its values by some constant.
Using the ranks (rather than constant values of $1,$ which is another option) weights the highest values where we want to focus attention.  This algorithm creates a running sum that goes up when a value from $A$ appears and (due to the negation) goes down when a value from $B$ appears.  If there's no real difference in their tails, this random walk should bounce up and down around zero.  (This is a consequence of the weighting by $1/m$ and $1/n.$)  If one of the tails is heavier, the random walk should initially trend upwards for a heavier $A$ tail and otherwise head downwards for a heavier $B$ tail.
This provides a nice diagnostic plot.  In the figure I have normalized the cumulative sum by multiplying all values by $1/\sqrt{n+m+1}$ and indexing them by the numbers $q = 0/(m+n), 1/(m+n), \ldots, (m+n-1)/(m+n).$ I call this the "cranksum" (cumulative rank sum).  Here is the first half, corresponding to the upper half of all the data:

There is a clear upward trend, consistent with what we saw in the previous figure.  But is it significant?
A simulation of the cranksums under the null hypothesis (of equally heavy tails) will settle this question.  Such a simulation creates many datasets of the same sizes as the original $A$ and $B$ (or, almost equivalently, creates many arbitrary permutations of the combined dataset) according to the same distribution (which distribution it is doesn't matter, provided it is continuous); computes their cranksums; and plots them.  Here are the first thousand out of 40,000 that I made for datasets of size $50$ and $100:$

The faint gray jagged curves in the middle form the assemblage of a thousand cranksum plots.  The yellow area, bounded in bold curves (the "envelope"), outlines the upper $99.25$ and lower $0.75$ percentiles of all 40,000 values.  Why these percentiles?  Because some analysis of these simulated data showed that only 5% of the simulated curves ever, at some point, go past these boundaries.  Thus, because the cranksum plot for the actual data does exceed the upper boundary for some of the initial (low) values of $q,$ it constitutes significant evidence at the $\alpha=0.05$ level that (1) the tails differ and (2) the tail of $A$ is heavier than the tail of $B.$
Of course you can see much more in the plot: the cranksum for our data is extremely high for all values of $q$ between $0$ and $0.23,$ approximately, and only then starts dropping, eventually reaching a height of $0$ around $q=0.5.$  Thus it is apparent that at least the upper $23\%$ of the underlying distribution of data set $A$ consistently exceeds the upper $23\%$ of the underlying distribution for dataset $B$ and likely the upper $50\%$ of ... $A$ exceeds the upper $50\%$ of ... $B.$
(Because these are synthetic data, I know their underlying distributions, so I can compute that for this example the CDFs cross at $x=1.2149$ at a height of $0.6515,$ implying the upper $34.85\%$ of the distribution for $A$ exceeds that of $B,$ quite in line with what the cranksum analysis is telling us based on the samples.)
Evidently it takes a little work to compute the cranksum and run the simulation, but it can be done efficiently: this simulation took two seconds, for instance.  To get you started, I have appended the R code used to make the figures.
#
# Testing whether one tail is longer than another.
# The return value is the cranksum, a vector of length m+n.
#
cranksum <- function(x, y) {
  m <- length(x)
  n <- length(y)
  i <- order(c(x,y))
  scores <- c(rep(1/m, m), rep(-1/n, n)) * rank(c(x,y))
  cumsum(scores[rev(i)]) / sqrt(n + m + 1)
}
#
# Create two datasets from two different distributions with the same means.
#
mu <- 0          # Logmean of `x`
sigma <- 1/2     # Log sd of `x`
k <- 20          # Gamma parameter of `y`
set.seed(17)
y <- rgamma(100, k, k/exp(mu + sigma^2/2)) # Gamma data
x <- exp(rnorm(50, mu, sigma))             # Lognormal data.
#
# Plot their ECDFs.
#
plot(ecdf(c(x,y)), cex=0, col="00000000", main="Empirical CDFs")
e.x <- ecdf(x)
curve(e.x(x), add=TRUE, lwd=2, n=1001)
e.y <- ecdf(y)
curve(e.y(x), add=TRUE, col="Red", lwd=2, n=1001)
#
# Simulate the null distribution (assuming no ties).
# Each simulated cranksum is in a column.
#
system.time(sim <- replicate(4e4, cranksum(runif(length(x)), runif(length(y)))))
#
# This alpha was found by trial and error, but that needs to be done only 
# once for any given pair of dataset sizes.
#
alpha <- 0.0075
tl <- apply(sim, 1, quantile, probs=c(alpha/2, 1-alpha/2)) # Cranksum envelope
#
# Compute the chances of exceeding the upper envelope or falling beneath the lower.
#
p.upper <- mean(apply(sim > tl[2,], 2, max))
p.lower <- mean(apply(sim < tl[1,], 2, max))
#
# Include the data with the simulation for the purpose of plotting everything together.
#
sim <- cbind(cranksum(x, y), sim)
#
# Plot.
#
q <- seq(0, 1, length.out=dim(sim)[1])
# The plot region:
plot(0:1/2, range(sim), type="n", xlab = "q", ylab = "Value", main="Cranksum Plot")
# The region between the envelopes:
polygon(c(q, rev(q)), c(tl[1,], rev(tl[2,])), border="Black", lwd=2, col="#f8f8e8")
# The cranksum curves themselves:
invisible(apply(sim[, seq.int(min(dim(sim)[2], 1e3))], 2, 
          function(y) lines(q, y, col="#00000004")))
# The cranksum for the data:
lines(q, sim[,1], col="#e01010", lwd=2)
# A reference axis at y=0:
abline(h=0, col="White")

A: I would suggest to fit different distributions on your observations, and to perform model selection to find the distribution that fits your observations the best. Exponential and Pareto distributions seem to be the best candidates given your hypotheses (positivity, monotone decrease). Once you have fitted these candidates distributions, model selection criteria, such as the Akaike Information Criterion (AIC) or the Bayesian Information Criterion (BIC) will give you a quantitative score for each model. The following paper will propose rules on how to interpret the evidences for the different models based on the BIC :
Kass, Robert E., and Adrian E. Raftery. "Bayes factors." Journal of the american statistical association 90.430 (1995): 773-795.
You may also want to have a look at this paper which deals with model inference of long-tailed distributions : 
Okada, Makoto, Kenji Yamanishi, and Naoki Masuda. "Long-tailed distributions of inter-event times as mixtures of exponential distributions." arXiv preprint arXiv:1905.00699 (2019). 
A: The OP wants a metric for "tail length." While that term is not precisely defined, one might assume that "tail heaviness" is desired. Both Pearson-based and quantile-based kurtosis are measures of tail heaviness.  (This useful application of kurtosis has long been unused because of the incorrect notion that kurtosis measures "peakedness" rather than "tail heaviness.")  See here https://math.stackexchange.com/questions/3521769/graphic-representation-of-kurtosis-and-skewness/3532888#3532888 for a clear explanation of why Pearson kurtosis measures tail heaviness.
Estimate such a tail heaviness by using the data for each sample, and find the sampling distribution of the difference.  Use this sampling distribution to assess (or test, as the OP wants) the "true magnitude" of the difference between the heaviness of the tails. You could use the bootstrap here, although some kinds of parametric or smoothed bootstrap analysis may be more reliable. (Tail heaviness is very difficult to estimate because it is only the rare extreme values in the data (or outliers) that provide the relevant information, and there are by definition very few of such data points.) 
