What is the best way to produce the actual pairs of points shown by QuantilePlot
? Would anyone have an implementation? (I am interested in the simplest possible version: QuantilePlot[l1, l2]
, where l1, l2
are lists of numbers.
2 Answers
Reverse-Engineering Mathematica's Q-Q-Plots
You can also use Quantile
directly which gives more control about what is happening and precision maybe:
Inspecting points produced by Bob Hanlon's approach and playing around with other sample sizes reveals that Mathematica appears to use the range specified by $\frac{k-0.3}{n+0.4}$ with $k = 1, \ldots, n$ to produce quantiles. Other possibilities are given by Wikipedia (cf. source no. 10 for the option given here).
Thus:
SeedRandom[0];
data1 = RandomVariate[NormalDistribution[2, 3], 100];
data2 = RandomVariate[StudentTDistribution[4, 2, 3], 200];
pts = With[
{
n = Min @@ Length /@ { data1, data2 }
},
Curry[Quantile][ Table[(k - 0.3)/(n + 0.4), {k, n}] ] /@ {data2, data1} // Transpose
(* or Map[ Quantile[#, Table[ ... ]]&] @ {data2,data1} // Transpose *)
];
ListPlot[ pts, Axes-> False, Frame -> True ]
We can compare with the QuantilePlot
points using Bob Hanlon's approach:
pts2 = Cases[plot, Point[pts_] :> pts, Infinity][[1]];
pts2 - pts
{{0.,0.}, ... , {0.,0.}}
This comparison also holds for samples sizes 200, 300 that I tested so far.
Thanks, Sjoerd C. de Vries, for pointing out the more general range-regime for Q-Q-Plots.
-
1$\begingroup$ Nice use of the new 11.3 function
Curry
. TheRange
call, though, is only correct because data 1 happens to be 100 points long. It doesn't generalize to differently sized data sets. $\endgroup$ Commented Mar 15, 2018 at 22:16 -
2$\begingroup$ Ah, you beat me to it. To add: the positioning method is due to Bernard and Bos-Levenbach (English version). Spelunking confirms that this method is used internally by the helper function
System`QuantilePlotDump`plottingPosition[]
. $\endgroup$ Commented Mar 16, 2018 at 12:15 -
1$\begingroup$ In the case where one is comparing against a distribution instead of another dataset,
EstimatedDistribution[]
is used on the data (NormalDistribution[]
is the default distribution to compare against), and then this is used in the computation of the second set of quantiles. $\endgroup$ Commented Mar 16, 2018 at 12:19
SeedRandom[0];
data1 = RandomVariate[NormalDistribution[2, 3], 100];
data2 = RandomVariate[StudentTDistribution[4, 2, 3], 200];
plt = QuantilePlot[data1, data2]
Use Cases
to extract the points
pts = Cases[plt, Point[pts_] :> pts, Infinity][[1]];
Dimensions[pts]
{100, 2}
In this case there are 100 points.
ListPlot[pts, Frame -> True, Axes -> False]
foo[[1, 2, 2, 2, 3, 1]]
works, but that is not the real question, since I have no idea what parameters mathematica uses (number of points, etc). I can guess, but this is not very satisfying. $\endgroup$QuantilePlot[]
. $\endgroup$