Memasang kurva kerapatan ke histogram di R

91

Apakah ada fungsi di R yang menyesuaikan kurva dengan histogram?

Katakanlah Anda memiliki histogram berikut

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

Kelihatannya normal, tapi miring. Saya ingin menyesuaikan kurva normal yang miring untuk membungkus histogram ini.

Pertanyaan ini agak mendasar, tetapi sepertinya saya tidak dapat menemukan jawaban untuk R di internet.

r histogram curve-fitting r-faq pengguna5243421
sumber
Apakah Anda ingin mencari m dan s sehingga distribusi Gaussian N (m, s) cocok dengan data Anda?
SteinNorheim
Saya tidak yakin apa artinya ...> _>
user5243421
10
@ mathee: Saya pikir maksudnya m = mean, dan s = deviasi standar. Distribusi gaussian adalah nama lain untuk distribusi normal.
Peter Mortensen

Jawaban:

154

Jika saya memahami pertanyaan Anda dengan benar, maka Anda mungkin menginginkan perkiraan kepadatan bersama dengan histogram:

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

Edit lama nanti:

Ini adalah versi yang sedikit lebih rapi:

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)

bersama dengan grafik yang dihasilkannya:

masukkan deskripsi gambar di sini

Dirk Eddelbuettel
sumber
3
+1 - dapatkah Anda juga melakukannya sebaliknya, yaitu menyesuaikan plot kepadatan agar sesuai dengan histogram?
vonjd
2
Saya sarankan memberikan parameter tambahan lines(density(X,na.rm= TRUE)karena vektor mungkin berisi nilai NA.
Anirudh
30

Hal seperti itu mudah dilakukan dengan 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()

atau untuk meniru hasil dari solusi Dirk

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

Inilah cara saya melakukannya:

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

Latihan bonus adalah melakukan ini dengan paket ggplot2 ...

John Johnson
sumber
Namun, jika Anda menginginkan sesuatu yang miring, Anda dapat melakukan contoh kepadatan dari atas, mengubah data Anda (mis. Foo.log & lt; - log (foo) dan coba yang di atas), atau coba paskan distribusi miring, seperti gamma atau lognormal (lognormal sama dengan mengambil log dan menyesuaikan dengan normal, btw).
John Johnson
2
Tetapi itu masih membutuhkan estimasi parameter distribusi Anda terlebih dahulu.
Dirk Eddelbuettel
Ini agak jauh dari sekadar membahas R, karena kita semakin memahami statistik teoretis, tetapi Anda dapat mencoba tautan ini untuk Gamma: en.wikipedia.org/wiki/Gamma_distribution#Parameter_estimation Untuk lognormal, ambil saja lognya (dengan asumsi semua data positif) dan bekerja dengan data yang ditransformasi log. Untuk apa pun yang lebih menarik, saya pikir Anda harus bekerja dengan buku teks statistik.
John Johnson
3
Saya pikir Anda salah paham bagaimana poster asli serta semua jawaban lainnya cukup puas menggunakan perkiraan non-parametrik - seperti histogram jadul atau perkiraan kepadatan berbasis data yang lebih modern. Perkiraan parametrik bagus jika Anda memiliki alasan kuat untuk mencurigai distribusi. Tapi bukan itu masalahnya di sini.
Dirk Eddelbuettel
11

Dirk telah menjelaskan cara memplot fungsi kerapatan di atas histogram. Namun terkadang Anda mungkin ingin menggunakan asumsi yang lebih kuat tentang distribusi normal miring dan plot yang bukannya kepadatan. Anda dapat memperkirakan parameter distribusi dan memplotnya menggunakan paket 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

Plot data terdistribusi miring-normal

Ini mungkin berfungsi lebih baik pada data yang lebih condong-normal:

Plot normal miring lainnya

fmark
sumber
3

Saya memiliki masalah yang sama tetapi solusi Dirk tampaknya tidak berhasil. Saya mendapatkan pesan peringatan ini setiap saat

"prob" is not a graphical parameter

Saya membaca ?histdan menemukan tentangfreq: a logical vector set TRUE by default.

kode yang berhasil untuk saya adalah

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