[R] ggplot2 의 커널 밀도 곡선(Kernel Density Curve)의 최대 피크값 좌표구하고 수직선 추가하기

이번 포스팅에서는 

(1) R ggplot2를 이용하여 커널 밀도 곡선(Kernel Density Curve)을 그리기, 

(2) 커밀 밀도 곡선의 최대 피크값의 좌표 구하기 
   (X, Y coordinates of the peak value in kernel density curve)

(3) 커널 밀도 곡선의 최대 피크값 위치에 수직선을 추가하는 방법 (adding a vertical line at peak value in kernel density plot) 을 소개하겠습니다. 

(1) R ggplot2로 커널 밀도 곡선 그리기 (Kernel Density Curve by R ggplot2)

예제로 사용할 데이터는 MASS 패키지에 내장되어 있는 Cars93 데이터프레임입니다. 자동차의 가격(Price 변수)에 대해서 R ggplot2로 커널 밀도 곡선 (not frequency, but density) 을 그려보겠습니다. 

> library(MASS) > library(ggplot2) >  > ggplot(Cars93, aes(x=Price)) +  +   geom_density(fill = 'yellow') +  +   geom_line(stat = "density") +  +   expand_limits(y = 0) +  +   ggtitle("Kernel Density Curve")

(2) 커널 밀도 곡선의 최대 피크값의 좌표 구하기       (X, Y coordinates of the peak value in kernel density curve)

density(df$var) 함수를 사용하면 x 와 y 값의 요약통계량 (최소, Q1, 중앙값, Q3, 최대값)을 확인할 수 있습니다. y 값의 최대값(Max)은 5.025e-02 값이네요. 

> density(Cars93$Price) Call: density.default(x = Cars93$Price) Data: Cars93$Price (93 obs.); Bandwidth 'bw' = 3.011        x                y              Min.   :-1.634   Min.   :1.608e-05    1st Qu.:16.508   1st Qu.:9.773e-04    Median :34.650   Median :5.329e-03    Mean   :34.650   Mean   :1.377e-02    3rd Qu.:52.792   3rd Qu.:2.078e-02    Max.   :70.934   Max.   :5.025e-02

이중에서 y 값 peak 값의 x 좌표를 알고 싶으므로 (a) which.max(density(Cars93$Price)$y) 로 y값의 최대 peak 값을 가지는 데이터의 index를 찾고 (127번째 데이터), 이 index를 이용해서 x 값의 좌표(16.25942)를 indexing 해올 수 있습니다. 

> which.max(density(Cars93$Price)$y) [1] 127 > x_coord_of_y_max = density(Cars93$Price)$x[127] > x_coord_of_y_max [1] 16.25942

즉, y (Price) 최대 피크값의 좌표는 (x, y) 는 (16.25922, 0.05025) 가 되겠습니다. 

(3) 커널 밀도 곡선의 최대 피크값 위치에 수직선을 추가하기      (adding a vertical line at peak value in kernel density plot)

geom_vline() 으로 (2)번에서 구한 y 피크값의 x 좌표를 입력해주면 커널밀도곡선의 최대 피크값 위치에 수직선(vertical line)을 추가할 수 있습니다. 

> ggplot(Cars93, aes(x=Price)) +  +   geom_density(fill = 'yellow') +  +   geom_line(stat = "density") +  +   expand_limits(y = 0) +  +   ggtitle("Kernel Density Curve w/ Vertical Line at Peak Point") + +   geom_vline(xintercept = x_coord_of_y_max, color='blue', size=2, linetype="dotted")

---

(4) 여러개 그룹의 커널 밀도 곡선을 겹쳐 그린 경우 최대 peak 값 찾고 수직선 그리기

위의 예는 1개의 데이터셋에 대해 1개의 커널 밀도 곡선을 그린 후 최대 peak 값을 찾는 것이었습니다. 이제부터는 여러개의 하위 그룹으로 나뉘어진 데이터셋으로 여러개의 커널 밀도 곡선을 겹쳐서 그린 경우에 최대 peak 값을 찾고, 그 위치에 빨간색으로 수직선을 추가해보겠습니다. 

예제로 Cars93 데이터프레임에서 차종(Type) 별 가격(Price)의 커널 밀도 곡선을 겹쳐서 그려보겠습니다. 

> ggplot(Cars93, aes(x=Price, colour = Type)) +  +   geom_density(fill = NA) +  +   geom_line(stat = "density") +  +   expand_limits(y = 0) +  +   ggtitle("Kernel Density Curve by Car Types")

차종(Type) 중에서 Van 의 커널 밀도 곡선이 최대 Peak 값에 해당하므로 --> 차종 중에서 Van 의 데이터만 가져와서 y 최대값과 이에 해당하는 관측치의 index 위치를 찾아보겠습니다. 

> density(Cars93[Cars93$Type=='Van', ]$Price) Call: density.default(x = Cars93[Cars93$Type == "Van", ]$Price) Data: Cars93[Cars93$Type == "Van", ]$Price (9 obs.); Bandwidth 'bw' = 0.303        x               y              Min.   :15.39   Min.   :0.0000072    1st Qu.:17.45   1st Qu.:0.0040424    Median :19.50   Median :0.0495649    Mean   :19.50   Mean   :0.1215343    3rd Qu.:21.55   3rd Qu.:0.1807290    Max.   :23.61   Max.   :0.5321582   > which.max(density(Cars93[Cars93$Type=='Van', 'Price'])$y) [1] 238 > x_coord_of_y_max = density(Cars93[Cars93$Type=='Van','Price'])$x[238] > x_coord_of_y_max [1] 19.20249

밀도 y의 최대값(max)은 0.532 이며, 이때의 관측치의 index는 238번째 값 이고, 관측치 238번째 값의 x값의 좌표는 19.20 입니다. 

따라서, 밀도 y 의 최대 peak 값의 좌표는 (x, y) = (19.20, 0.53) 이네요. 

마지막으로, geom_vline() 을 사용해서 밀도 최대 peak 값에 해당하는 x=19.20 을 기준으로 빨간색 점선 수직선(vertical line)을 추가해보겠습니다. 

> ggplot(Cars93, aes(x=Price, colour = Type)) +  +   geom_density(fill = NA) +  +   geom_line(stat = "density") +  +   expand_limits(y = 0) +  +   ggtitle("Kernel Density Curve w/ Vertical Line at Peak Point") + +   theme_bw() + +   geom_vline(xintercept = x_coord_of_y_max, color='red', size=1, linetype="dashed")

R ggplot2로 히스토그램, 커널밀도곡선을 그리는 다양한 방법은 https://rfriend.tistory.com/67 를 참고하세요. 

많은 도움이 되었기를 바랍니다. 

이번 포스팅이 도움이 되었다면 아래의 '공감~'를 꾹 눌러주세요. ^^