Comment interpréter un graphique QQ

Je travaille avec un petit ensemble de données (21 observations) et ai le graphe QQ normal suivant dans R:

Voyant que l'intrigue ne soutient pas la normalité, que pourrais-je déduire de la distribution sous-jacente? Il me semble qu'une distribution plus biaisée à droite conviendrait mieux, n'est-ce pas? Aussi, quelles autres conclusions pouvons-nous tirer des données?

Si les valeurs sont situées le long d'une ligne, la distribution a la même forme (jusqu'à l'emplacement et l'échelle) que la distribution théorique que nous avons supposée.

Comportement local : Lorsque vous examinez les valeurs d’échantillon triées sur l’axe des y et les quantiles attendus (approximatifs) sur l’axe des x, vous pouvez identifier la différence entre les valeurs d’une section du tracé et une tendance linéaire globale en vérifiant si le les valeurs sont plus ou moins concentrées que ne le supposerait la distribution théorique dans cette partie de la parcelle:

Comme nous le voyons, les points moins concentrés augmentent de plus en plus de points que supposés augmenter moins rapidement que ne le suggère une relation linéaire globale et, dans les cas extrêmes, correspondent à un écart dans la densité de l'échantillon (se présente sous la forme d'un saut presque vertical) ou un pic de valeurs constantes (valeurs alignées horizontalement). Cela nous permet de repérer une queue épaisse ou légère et, partant, une asymétrie supérieure ou inférieure à la distribution théorique, etc.

Apparence générale:

Voici à quoi ressemblent les tracés QQ (pour des choix particuliers de distribution) en moyenne :

Mais le hasard a tendance à obscurcir les choses, surtout avec de petits échantillons:

Notez que, à les résultats peuvent être beaucoup plus variables que ceux présentés ici. J'ai généré plusieurs ensembles de six parcelles de ce type et choisi un ensemble "sympa" dans lequel vous pouvez voir la forme des six parcelles en même temps. Parfois, les relations droites semblent incurvées, les relations incurvées semblent rectilignes, les queues épaisses semblent obliques, et ainsi de suite - avec de si petits échantillons, la situation peut souvent être beaucoup moins claire:$n=21$

Il est possible de distinguer plus de caractéristiques que celles (comme la discrétion, par exemple), mais avec , même de telles caractéristiques de base peuvent être difficiles à repérer; nous ne devrions pas essayer de «surinterpréter» chaque petit mouvement. À mesure que la taille des échantillons augmente, les parcelles se «stabilisent» en général et les caractéristiques deviennent plus facilement interprétables que représentatives du bruit. [Avec certaines distributions très lourdes, la grande valeur rare rare pourrait empêcher la bonne stabilisation de l'image, même avec des échantillons de taille assez grande.]$n=21$

Vous pouvez également trouver la suggestion ici utile lorsque vous essayez de décider de combien vous devriez vous inquiéter d’une quantité particulière de courbure ou de wiggliness.

Un guide plus approprié pour l'interprétation en général comprendrait également des affichages pour des échantillons de plus en plus grands.

J'ai fait une application brillante pour aider à interpréter l'intrigue QQ normale. Essayez ce lien.

Dans cette application, vous pouvez régler l'asymétrie, la finesse (kurtosis) et la modalité de données et vous pouvez voir comment l'histogramme et le tracé de QQ changent. Inversement, vous pouvez l'utiliser d'une manière qui, compte tenu du schéma du tracé QQ, vérifie ensuite comment doit être l'asymétrie, etc.

Pour plus de détails, voir la documentation correspondante.

J'ai réalisé que je n'avais pas assez d'espace libre pour fournir cette application en ligne. Comme la demande, je fournirai tous les trois morceaux de code: sample.R, server.Ret ui.Rici. Ceux qui sont intéressés par l'exécution de cette application peuvent simplement charger ces fichiers dans Rstudio puis l'exécuter sur votre propre PC.

Le sample.Rfichier:

# Compute the positive part of a real number x, which is $\max(x, 0)$.
positive_part <- function(x) {ifelse(x > 0, x, 0)}

# This function generates n data points from some unimodal population.
# Input: ----------------------------------------------------
# n: sample size;
# mu: the mode of the population, default value is 0.
# skewness: the parameter that reflects the skewness of the distribution, note it is not
#           the exact skewness defined in statistics textbook, the default value is 0.
# tailedness: the parameter that reflects the tailedness of the distribution, note it is
#             not the exact kurtosis defined in textbook, the default value is 0.

# When all arguments take their default values, the data will be generated from standard 
# normal distribution.

random_sample <- function(n, mu = 0, skewness = 0, tailedness = 0){
  sigma = 1

  # The sampling scheme resembles the rejection sampling. For each step, an initial data point
  # was proposed, and it will be rejected or accepted based on the weights determined by the
  # skewness and tailedness of input. 
  reject_skewness <- function(x){
      scale = 1
      # if `skewness` > 0 (means data are right-skewed), then small values of x will be rejected
      # with higher probability.
      l <- exp(-scale * skewness * x)
      l/(1 + l)
  }

  reject_tailedness <- function(x){
      scale = 1
      # if `tailedness` < 0 (means data are lightly-tailed), then big values of x will be rejected with
      # higher probability.
      l <- exp(-scale * tailedness * abs(x))
      l/(1 + l)
  }

  # w is another layer option to control the tailedness, the higher the w is, the data will be
  # more heavily-tailed. 
  w = positive_part((1 - exp(-0.5 * tailedness)))/(1 + exp(-0.5 * tailedness))

  filter <- function(x){
    # The proposed data points will be accepted only if it satified the following condition, 
    # in which way we controlled the skewness and tailedness of data. (For example, the 
    # proposed data point will be rejected more frequently if it has higher skewness or
    # tailedness.)
    accept <- runif(length(x)) > reject_tailedness(x) * reject_skewness(x)
    x[accept]
  }

  result <- filter(mu + sigma * ((1 - w) * rnorm(n) + w * rt(n, 5)))
  # Keep generating data points until the length of data vector reaches n.
  while (length(result) < n) {
    result <- c(result, filter(mu + sigma * ((1 - w) * rnorm(n) + w * rt(n, 5))))
  }
  result[1:n]
}

multimodal <- function(n, Mu, skewness = 0, tailedness = 0) {
  # Deal with the bimodal case.
  mumu <- as.numeric(Mu %*% rmultinom(n, 1, rep(1, length(Mu))))
  mumu + random_sample(n, skewness = skewness, tailedness = tailedness)
}

Le server.Rfichier:

library(shiny)
# Need 'ggplot2' package to get a better aesthetic effect.
library(ggplot2)

# The 'sample.R' source code is used to generate data to be plotted, based on the input skewness, 
# tailedness and modality. For more information, see the source code in 'sample.R' code.
source("sample.R")

shinyServer(function(input, output) {
  # We generate 10000 data points from the distribution which reflects the specification of skewness,
  # tailedness and modality. 
  n = 10000

  # 'scale' is a parameter that controls the skewness and tailedness.
  scale = 1000

  # The `reactive` function is a trick to accelerate the app, which enables us only generate the data
  # once to plot two plots. The generated sample was stored in the `data` object to be called later.
  data <- reactive({
    # For `Unimodal` choice, we fix the mode at 0.
    if (input$modality == "Unimodal") {mu = 0}

    # For `Bimodal` choice, we fix the two modes at -2 and 2.
    if (input$modality == "Bimodal") {mu = c(-2, 2)}

    # Details will be explained in `sample.R` file.
    sample1 <- multimodal(n, mu, skewness = scale * input$skewness, tailedness = scale * input$kurtosis)
    data.frame(x = sample1)})

  output$histogram <- renderPlot({
    # Plot the histogram.
    ggplot(data(), aes(x = x)) + 
      geom_histogram(aes(y = ..density..), binwidth = .5, colour = "black", fill = "white") + 
      xlim(-6, 6) +
      # Overlay the density curve.
      geom_density(alpha = .5, fill = "blue") + ggtitle("Histogram of Data") + 
      theme(plot.title = element_text(lineheight = .8, face = "bold"))
  })

  output$qqplot <- renderPlot({
    # Plot the QQ plot.
    ggplot(data(), aes(sample = x)) + stat_qq() + ggtitle("QQplot of Data") + 
      theme(plot.title = element_text(lineheight=.8, face = "bold"))
    })
})

Enfin, le ui.Rfichier:

library(shiny)

# Define UI for application that helps students interpret the pattern of (normal) QQ plots. 
# By using this app, we can show students the different patterns of QQ plots (and the histograms,
# for completeness) for different type of data distributions. For example, left skewed heavy tailed
# data, etc. 

# This app can be (and is encouraged to be) used in a reversed way, namely, show the QQ plot to the 
# students first, then tell them based on the pattern of the QQ plot, the data is right skewed, bimodal,
# heavy-tailed, etc.


shinyUI(fluidPage(
  # Application title
  titlePanel("Interpreting Normal QQ Plots"),

  sidebarLayout(
    sidebarPanel(
      # The first slider can control the skewness of input data. "-1" indicates the most left-skewed 
      # case while "1" indicates the most right-skewed case.
      sliderInput("skewness", "Skewness", min = -1, max = 1, value = 0, step = 0.1, ticks = FALSE),

      # The second slider can control the skewness of input data. "-1" indicates the most light tail
      # case while "1" indicates the most heavy tail case.
      sliderInput("kurtosis", "Tailedness", min = -1, max = 1, value = 0, step = 0.1, ticks = FALSE),

      # This selectbox allows user to choose the number of modes of data, two options are provided:
      # "Unimodal" and "Bimodal".
      selectInput("modality", label = "Modality", 
                  choices = c("Unimodal" = "Unimodal", "Bimodal" = "Bimodal"),
                  selected = "Unimodal"),
      br(),
      # The following helper information will be shown on the user interface to give necessary
      # information to help users understand sliders.
      helpText(p("The skewness of data is controlled by moving the", strong("Skewness"), "slider,", 
               "the left side means left skewed while the right side means right skewed."), 
               p("The tailedness of data is controlled by moving the", strong("Tailedness"), "slider,", 
                 "the left side means light tailed while the right side means heavy tailedd."),
               p("The modality of data is controlledy by selecting the modality from", strong("Modality"),
                 "select box.")
               )
  ),

  # The main panel outputs two plots. One plot is the histogram of data (with the nonparamteric density
  # curve overlaid), to get a better visualization, we restricted the range of x-axis to -6 to 6 so 
  # that part of the data will not be shown when heavy-tailed input is chosen. The other plot is the 
  # QQ plot of data, as convention, the x-axis is the theoretical quantiles for standard normal distri-
  # bution and the y-axis is the sample quantiles of data. 
  mainPanel(
    plotOutput("histogram"),
    plotOutput("qqplot")
  )
)
)
)

Une explication très utile (et intuitive) est donnée par prof. Philippe Rigollet dans le cours MOOC MIT: 18.650 Statistiques pour les applications, automne 2016 - voir la vidéo à 45 minutes

https://www.youtube.com/watch?v=vMaKx9fmJHE

J'ai copié grossièrement son diagramme que je garde dans mes notes car je le trouve très utile.

Dans l'exemple 1, dans le diagramme en haut à gauche, nous voyons que dans la queue droite le quantile empirique (ou échantillon) est inférieur au quantile théorique

Qe <Qt

Ceci peut être interprété à l'aide des fonctions de densité de probabilité. Pour la même valeur , le quantile empirique se trouve à gauche du quantile théorique, ce qui signifie que la queue droite de la distribution empirique est "plus légère" que la queue droite de la distribution théorique, c’est-à-dire qu’elle tombe plus rapidement à des valeurs proches de zéro.$\alpha$

Étant donné que ce fil a été jugé définitif par StackExchange "comment interpréter le tracé de qq normal", je voudrais indiquer aux lecteurs une relation mathématique précise et agréable entre le tracé de qq normal et la statistique de kurtosis en excès.

C'est ici:

https://stats.stackexchange.com/a/354076/102879

Un résumé bref (et trop simplifié) est donné comme suit (voir le lien pour des énoncés mathématiques plus précis): Vous pouvez réellement voir un kurtosis en excès dans le tracé normal de qq comme la distance moyenne entre les quantiles de données et les quantiles normaux théoriques correspondants, pondérée. par distance des données à la moyenne. Ainsi, lorsque les valeurs absolues dans les queues du graphique qq s'écartent généralement des valeurs normales prévues dans les directions extrêmes, vous obtenez un excès de kurtosis.

Comme kurtosis est la moyenne de ces déviations pondérée par les distances par rapport à la moyenne, les valeurs proches du centre du graphique qq ont peu d’impact sur kurtosis. Par conséquent, l'excès de kurtosis n'est pas lié au centre de la distribution, où se trouve le "pic". L'excès de kurtosis est plutôt presque entièrement déterminé par la comparaison des queues de la distribution de données à la distribution normale.