Ajustement d'une courbe de densité à un histogramme dans R

91

Y a-t-il une fonction dans R qui adapte une courbe à un histogramme?

Disons que vous avez l'histogramme suivant

hist(c(rep(65, times=5), rep(25, times=5), rep(35, times=10), rep(45, times=4)))

Cela semble normal, mais il est biaisé. Je veux ajuster une courbe normale qui est biaisée pour s'enrouler autour de cet histogramme.

Cette question est assez basique, mais je n'arrive pas à trouver la réponse pour R sur Internet.

r histogram curve-fitting r-faq user5243421
la source
Voulez-vous trouver m et s tels que la distribution gaussienne N (m, s) s'adapte à vos données?
SteinNorheim
Je ne sais pas ce que cela signifie ...> _>
user5243421
10
@mathee: Je pense qu'il veut dire m = moyenne et s = écart type. La distribution gaussienne est un autre nom pour la distribution normale.
Peter Mortensen

Réponses:

154

Si je comprends bien votre question, vous voulez probablement une estimation de la densité avec l'histogramme:

X <- c(rep(65, times=5), rep(25, times=5), rep(35, times=10), rep(45, times=4))
hist(X, prob=TRUE)            # prob=TRUE for probabilities not counts
lines(density(X))             # add a density estimate with defaults
lines(density(X, adjust=2), lty="dotted")   # add another "smoother" density

Modifier longtemps plus tard:

Voici une version un peu plus habillée:

X <- c(rep(65, times=5), rep(25, times=5), rep(35, times=10), rep(45, times=4))
hist(X, prob=TRUE, col="grey")# prob=TRUE for probabilities not counts
lines(density(X), col="blue", lwd=2) # add a density estimate with defaults
lines(density(X, adjust=2), lty="dotted", col="darkgreen", lwd=2)

avec le graphique qu'il produit:

entrez la description de l'image ici

Dirk Eddelbuettel
la source
3
+1 - pouvez-vous également le faire dans l'autre sens, c'est-à-dire ajuster le graphique de densité pour l'adapter à l'histogramme?
vonjd
2
Je suggère de donner un paramètre supplémentaire à lines(density(X,na.rm= TRUE)car le vecteur peut contenir des valeurs NA.
Anirudh le
30

Une telle chose est facile avec ggplot2

library(ggplot2)
dataset <- data.frame(X = c(rep(65, times=5), rep(25, times=5), 
                            rep(35, times=10), rep(45, times=4)))
ggplot(dataset, aes(x = X)) + 
  geom_histogram(aes(y = ..density..)) + 
  geom_density()

ou pour imiter le résultat de la solution de Dirk

ggplot(dataset, aes(x = X)) + 
  geom_histogram(aes(y = ..density..), binwidth = 5) + 
  geom_density()
Thierry
la source
28

Voici comment je le fais:

foo <- rnorm(100, mean=1, sd=2)
hist(foo, prob=TRUE)
curve(dnorm(x, mean=mean(foo), sd=sd(foo)), add=TRUE)

Un exercice bonus est de faire cela avec le package ggplot2 ...

John Johnson
la source
Cependant, si vous voulez quelque chose qui est biaisé, vous pouvez soit faire l'exemple de densité ci-dessus, transformer vos données (par exemple foo.log & lt; - log (foo) et essayer ce qui précède), soit essayer d'ajuster une distribution biaisée, telle que le gamma ou lognormal (lognormal équivaut à prendre le journal et à ajuster une normale, btw).
John Johnson
2
Mais cela nécessite toujours d'estimer d'abord les paramètres de votre distribution.
Dirk Eddelbuettel
Cela est un peu loin de la simple discussion de R, car nous abordons davantage les statistiques théoriques, mais vous pouvez essayer ce lien pour le Gamma: en.wikipedia.org/wiki/Gamma_distribution#Parameter_estimation Pour lognormal, prenez simplement le journal (en supposant toutes les données sont positives) et fonctionnent avec des données transformées en log. Pour tout ce qui est plus sophistiqué, je pense que vous devriez travailler avec un manuel de statistiques.
John Johnson
3
Je pense que vous comprenez mal comment l'affiche originale ainsi que toutes les autres réponses se contentent d'utiliser des estimations non paramétriques - comme un histogramme à l'ancienne ou une estimation de densité basée sur des données un peu plus moderne. Les estimations paramétriques sont excellentes si vous avez de bonnes raisons de soupçonner une distribution. Mais ce n'était pas le cas ici.
Dirk Eddelbuettel
11

Dirk a expliqué comment tracer la fonction de densité sur l'histogramme. Mais parfois, vous voudrez peut-être aller avec l'hypothèse plus forte d'une distribution normale asymétrique et tracer cela au lieu de la densité. Vous pouvez estimer les paramètres de la distribution et la tracer à l'aide du package sn :

> sn.mle(y=c(rep(65, times=5), rep(25, times=5), rep(35, times=10), rep(45, times=4)))
$call
sn.mle(y = c(rep(65, times = 5), rep(25, times = 5), rep(35, 
    times = 10), rep(45, times = 4)))

$cp
    mean     s.d. skewness 
41.46228 12.47892  0.99527

Diagramme de données distribuées asymétrique

Cela fonctionne probablement mieux sur des données plus faussées-normales:

Un autre tracé asymétrique

fmark
la source
3

J'ai eu le même problème mais la solution de Dirk ne semblait pas fonctionner. Je recevais ce message d'avertissement à chaque fois

"prob" is not a graphical parameter

J'ai lu ?histet trouvé surfreq: a logical vector set TRUE by default.

le code qui a fonctionné pour moi est

hist(x,freq=FALSE)
lines(density(x),na.rm=TRUE)
Matias Andina
la source