Best way to construct a QQ-plot

Question

I want to assess the normality of a dataset (which is log-normally distributed data transformed back to normal) using a Q-Q plot.

I stumbled on the fact that there are many ways to build such a plot, as there are multiple ways to determine sample quantiles; and different ways to select where to place them with respect to the theoretical quantiles.

As an illustration, I generated a random sample of size $N=50$ and plotted the quantiles of this sample against theoretical quantiles. The blue dots correspond to theoretical and sample quantiles $q_a = \{1,\ldots, N\}/N$; the orange dots correspond to theoretical quantiles as per Filliben's estimate, while the matching sample quantiles are still $q_a$. The figures' titles corresponds to the way the sample quantile have been determined; R* notation references this table; HD references the Harrell-Davis estimator, which is not mentioned on the Wikipedia page. SciPy default corresponds to the default behaviour of scipy.mstats.mquantiles, which does not seem to be documented on Wikipedia either.

The orange dots of the R4 plot correspond to the behaviour of scipy.stats.problot, which plots the sorted data against the theoretical quantiles evaluated with the Filliben's estimate.

This second figure shows these different quantile estimates against $q_a$.

The differences are rather marginal. The HD estimator seems to yield a smoother curve, which tempts me in using it, but it is a rather shallow justification.

I obtain my actual data by combining several datasets with different weights; I will definitely have more than 50 points, likely ranging from $\mathcal{O}(10^{4})$ to $\mathcal{O}(10^6)$. Different quantile estimators might lead to different methods to incorporate these weights, hence why I'd like to make an "educated" choice of said estimator before working out this problem. Any input is appreciated.

As a final precision: I have to assume the data is log-normally distributed, as I need to estimate the distribution's parameters to scale it back to a standard normal distribution before building the Q-Q plots. It's more an a posteriori verification of what had been hinted by histograms than actual exploration.

EDIT: I realized that I don't really need all the distribution parameters, only the shift before taking the log to compare with a standard normal. In my case it can reasonably be anything between 0 and the first order statistic.

---

---

To expand on my actual problem, here are some examples of my data. The "Sample quantiles" axis correspond to data transformed to standard normal, i.e. $X = \frac{\ln(Y - \tau) - \mu}{\sigma}$ where $Y$ would be the actual data. Data points are in blue, the orange line goes through the first and third quartiles and the black points are 20 realisations of a random variable sampled from the standard normal distribution, as suggested in @BruceET's answer.

The key here is my data points are weighted. I cannot just sort them to get sample quantiles; my quantiles are a linear combination of the data points, this transformation depending on the weights. The matrix of the transformation is quite sparse, but I have not yet managed to build it more efficiently than linear time with respect to the number of quantiles I am interested in. Therefore, instead of determining the quantiles against the theoretical ones evaluated at, say, seq(.5/n, 1-.5/n, length=n) with $n=N$ and $N$ my number of observations, I do so with a way smaller $n$.

This plot clearly suggests my data ($N \approx 4\cdot 10^6$ is log-normally distributed with correctly estimated parameters.

This one clearly suggests either the data ($N \approx 10^7$) is not log-normally distributed or the estimation has failed.

In both cases, I use $n=1001$.

The way I determined the "quantile weights" from my actual weights is based on the method described here, using the R7 definition. Note that the Filliben's estimate goes out of the window when not working with the sorted original dataset.

I now see more clearly that there is not a method for quantile estimation; however, I am wondering to what extent not using all the data points impacts the resulting plot.

Answer 1

Your normal Q-Q plots look fine. For a sample of size $n=50,$ it seems you may be striving for more precision than Q-Q plots usually provide.

Several styles of normal Q-Q plots are used to judge normality of a possibly normal sample. Some Q-Q plots put data on the horizontal axis and some on the vertical axis. (Preferred usage seems to vary by country.) In R, the default is to put data on the vertical axis, but the parameter datax=T will put data on the horizontal axis.

Looking at plots is sometimes more useful than formal tests of goodness-of-fit to a normal distribution. But looking at a plot is not the same thing as a formal test, so assessing normality from a plot requires familiarity with the kind of plot you use. The general idea is that the normal quantile function is transformed to a straight line in making a Q-Q plot. So the Q-Q plot of a normal sample should be 'nearly' linear (perhaps with more deviation from a line in the tails than in the center of the sample).

Two kinds of guide lines have been used with the Q-Q plots made in R. One kind of line (blue, left below) connects the lower quartiles (data and theoretical) with the upper quartiles. The other uses $y = -\bar X/S + x/S,$ where $\bar X$ and $S$ are the sample mean and SD, respectively (when data is on the horizontal axis).

set.seed(2021)
x = rnorm(100, 50, 7)  # SD is 7
par(mfrow=c(1,2))
 qqnorm(x, datax=T)
  qqline(x, datax=T, col="blue", lwd=2)
 qqnorm(x, datax=T)
  abline(-mean(x)/sd(x), 1/sd(x), col="red", lwd=2)
par(mfrow=c(1,1))

Some statistical programs put curved 'confidence bands' around the data cloud made by the Q-Q plot, with the idea that "too many" data points outside the bands may serve to reject in an informal test of normality.

Another approach is to superimpose 20 additional Q-Q plots in the background of the Q-Q plot of the data, to suggest an informal 95% confidence region. Data known to be sampled from a normal distribution are used. The mean and variance of these additional samples may be known or hypothetical $\mu$ and $\sigma$ or the mean and SD of the sample being plotted.

To make these additional Q-Q plots (and for other purposes) it is sometimes convenient to be able to make Q-Q plots 'from scratch' instead of using a pre-programmed procedure.

For a sample of size $n$ it is convenient to use theoretical quantiles and sorted data as below. The figure below uses the same data x as above.

th.quant = qnorm(seq(.5/n, 1-.5/n, length=n))
dta.quant = sort(x)

qqnorm(x, datax=T)
n = 100;  th.quant = qnorm(seq(.5/n, 1-.5/n, length=n))
for(i in 1:20) {
 points(sort(rnorm(n,mean(x),sd(x))), th.quant, col="green")
}
points(sort(x),th.quant, pch=19)  # refresh