'R 분석과 프로그래밍' 카테고리의 글 목록 (20 Page)

'R 분석과 프로그래밍'에 해당되는 글 289건

2015.10.15 R (1) 독립된 두 표본의 모평균 차이에 대한 추정과 검정 : t.test() 11
2015.10.12 R 두 연속형 변수 상관관계 : 공분산(covariance) cov(), 상관계수(correlation coefficient) cor() 4
2015.10.12 R 연속형 변수 (2) 그룹별(요인별) 요약통계 비교하기 : tapply(var, factor, summary), by(), aggregate(), summaryBy(), describeBy() 4
2015.10.12 R 연속형 변수 (1) 요약통계 한번에 보기 : summary(), stat.desc(), describe() 2
2015.10.12 R 일변량 연속형 자료 요약 - (3) 분포형태와 대칭정도(distrubution, symmetry) : 왜도(skewness), 첨도(kurtosis), 분위수-분위수(Quantile-Qunantile)
2015.10.11 R 일변량 연속형 자료 요약 - (2) 퍼짐 정도(dispersion) : 분산(variance), 표준편차(sd), 변이계수, 범위(range), IQR, 백분위수(percentile)
2015.10.11 R 일변량 연속형 자료 요약 - (1) 중심화 경향 : 산술평균(mean), 중앙값(median), 최빈값, 기하평균, CAGR, 조화평균, 가중평균 2
2015.10.10 R 이변량 범주형 자료의 정리와 요약 - 분할표(contingency table) : table(), xtabs(), prop.table(), margin.table(), addmargins(), CrossTable()
2015.10.10 R 일변량 범주형 자료의 정리와 요약 - 도수분포표(frequency distribution table) : table(), xtabs() 2
2015.10.09 R 단일 모집단 분포의 정규성 검정 : shapiro.test(), qqnorm(), qqline()

R (1) 독립된 두 표본의 모평균 차이에 대한 추정과 검정 : t.test()

R 분석과 프로그래밍/R 통계분석 2015. 10. 15. 01:05

두개의 모집단에 대한 추정과 검정 (two sample tests)에 대해서 정규분포 가정을 만족하는 경우와 그렇지 않은 경우로 나누어서 알아보겠습니다.

정규성 가정을 충족하는 경우

(1) 독립된 두 표본의 평균 차이에 대한 추정과 검정 : t.test()

(indepentent two sample t-test)

(2) 짝을 이룬 표본에 대한 평균 차이에 대한 추정과 검정 : t.test(paired=TRUE)

(paired sample t-test)

(3) 두 모집단의 모비율 차이에 대한 추정과 검정 : prop.test()

(independent two population proportions test)

정규성 가정을 충족하지 못하는 경우, 혹은 분포형태를 모르는 경우

(4) 두 모집단의 중심 차이에 대한 비모수 검정 : wilcox.test()

(non-parametric wilcoxon tests on two indepedent sample)

먼저, (1) 독립된 두 표본의 평균 차이에 대한 추정과 검정(indepentent two sample t-test)를 R의 t.test() 함수를 이용해 분석해 보겠습니다.

예제로 사용할 데이터는 MASS 패키지에 내장된 Cars93 데이터프레임의 가격(Price)과 생산국가(Origin) 입니다. 생산국이 USA vs. non-USA 2개의 group 에 대해서 차 가격(Price)의 평균이 차이가 있는지를 검정해보고, 95% 신뢰구간을 추정해보겠습니다.

> # Dataset
> library(MASS)
> str(Cars93)
'data.frame':	93 obs. of  27 variables:
 $ Manufacturer      : Factor w/ 32 levels "Acura","Audi",..: 1 1 2 2 3 4 4 4 4 5 ...
 $ Model             : Factor w/ 93 levels "100","190E","240",..: 49 56 9 1 6 24 54 74 73 35 ...
 $ Type              : Factor w/ 6 levels "Compact","Large",..: 4 3 1 3 3 3 2 2 3 2 ...
 $ Min.Price         : num  12.9 29.2 25.9 30.8 23.7 14.2 19.9 22.6 26.3 33 ...
 $ Price             : num  15.9 33.9 29.1 37.7 30 15.7 20.8 23.7 26.3 34.7 ...
 $ Max.Price         : num  18.8 38.7 32.3 44.6 36.2 17.3 21.7 24.9 26.3 36.3 ...
 $ MPG.city          : int  25 18 20 19 22 22 19 16 19 16 ...
 $ MPG.highway       : int  31 25 26 26 30 31 28 25 27 25 ...
 $ AirBags           : Factor w/ 3 levels "Driver & Passenger",..: 3 1 2 1 2 2 2 2 2 2 ...
 $ DriveTrain        : Factor w/ 3 levels "4WD","Front",..: 2 2 2 2 3 2 2 3 2 2 ...
 $ Cylinders         : Factor w/ 6 levels "3","4","5","6",..: 2 4 4 4 2 2 4 4 4 5 ...
 $ EngineSize        : num  1.8 3.2 2.8 2.8 3.5 2.2 3.8 5.7 3.8 4.9 ...
 $ Horsepower        : int  140 200 172 172 208 110 170 180 170 200 ...
 $ RPM               : int  6300 5500 5500 5500 5700 5200 4800 4000 4800 4100 ...
 $ Rev.per.mile      : int  2890 2335 2280 2535 2545 2565 1570 1320 1690 1510 ...
 $ Man.trans.avail   : Factor w/ 2 levels "No","Yes": 2 2 2 2 2 1 1 1 1 1 ...
 $ Fuel.tank.capacity: num  13.2 18 16.9 21.1 21.1 16.4 18 23 18.8 18 ...
 $ Passengers        : int  5 5 5 6 4 6 6 6 5 6 ...
 $ Length            : int  177 195 180 193 186 189 200 216 198 206 ...
 $ Wheelbase         : int  102 115 102 106 109 105 111 116 108 114 ...
 $ Width             : int  68 71 67 70 69 69 74 78 73 73 ...
 $ Turn.circle       : int  37 38 37 37 39 41 42 45 41 43 ...
 $ Rear.seat.room    : num  26.5 30 28 31 27 28 30.5 30.5 26.5 35 ...
 $ Luggage.room      : int  11 15 14 17 13 16 17 21 14 18 ...
 $ Weight            : int  2705 3560 3375 3405 3640 2880 3470 4105 3495 3620 ...
 $ Origin            : Factor w/ 2 levels "USA","non-USA": 2 2 2 2 2 1 1 1 1 1 ...
 $ Make              : Factor w/ 93 levels "Acura Integra",..: 1 2 4 3 5 6 7 9 8 10 ...

(1) 탐색적 분석

t-test에 들어가기 전에 탐색적 데이터 분석의 일환으로 R로 Origin별 관측값의 개수를 살펴보고, summary statistics, 박스그림(Box plot)과 Histogram 을 그려보면 아래와 같습니다.

- USA 48개, non-USA 45개 관측치로서 sample size가 작지는 않으므로 t-test에 문제는 없겠습니다.

- summary 통계량(mean), Box plot과 Histogram을 보니 두 집단(USA vs. non-USA)의 가격의 평균에는 큰 차이는 없어보이는군요.

- summary 통계량(Min, Max), Box plot과 Histogram을 보니 두 집단(USA vs. non-USA)의 가격의 분산은 아무래도 차이가 있어보이는데요,

이것을 t-검정(t-test)을 통해서 확인해보겠습니다.

> # frequency distribution table
> table(Cars93$Origin)

    USA non-USA 
     48      45 
> 
>

> # summary
> with(Cars93, tapply(Price, Origin, summary))
$USA
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
   7.40   13.48   16.30   18.57   20.72   40.10 

$`non-USA`
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
   8.00   11.60   19.10   20.51   26.70   61.90
> 
> 
> # box plot
> boxplot(Price ~ Origin, # Price by Origin (USA, non-USA)
+         data = Cars93, 
+         main = "Boxplot of Price by Origin",  
+         xlab = "Origin",  # x axis label
+         ylab = "Price") # y axis label
>

> 
> # Histogram
> # install.packages("ggplot2") # only for first user of ggplot2
> library(ggplot2)
> 
> ggplot(Cars93, aes(x=Price)) + 
+   geom_histogram(binwidth=5) + 
+   facet_grid(Origin ~ .) +
+   ggtitle("Histogram of Price by Origin")

(2) 두 모집단의 분산 동일성 가정에 대한 검정 : F-test

두 독립 표본의 t-test 에서 두 모집단의 분산 동일성 가정을 만족하는지 여부에 따라 분석 옵션을 달리 주어야 하므로, 먼저 분산 동일성에 대한 F-검정(F-test)을 var.test() 함수를 이용해 해보겠습니다.

> # independent two sample variance test : var.test()
> # H0 : variance of price by origin is equal (sigma1^2/sigma2^2  = 1)
> # H1 : variance of price by origin is not equal (sigma1^2/sigma2^2  != 1)
> 
> var.test(Price ~ Origin, data = Cars93)

	F test to compare two variances

data:  Price by Origin
F = 0.47796, num df = 47, denom df = 44, p-value = 0.01387
alternative hypothesis: true ratio of variances is not equal to 1
95 percent confidence interval:
 0.2645004 0.8587304
sample estimates:
ratio of variances 
         0.4779637

P-value가 0.01387로서 매우 작은 값이 나왔으므로, 유의수준 5% 하에서 두 모집단의 분산이 동일하다는 귀무가설을 기각하고, 분산이 서로 다르다는 대립가설을 채택하게 됩니다. 따라서 아래에 t-test 할 때는 분산 동일성 여부에 대한 옵션에서 FALSE (var.equal = FALSE)를 지정하면 되겠습니다.

R의 t-test() 함수를 이용하는 방법에는 2가지가 있습니다. 첫번째는 x (numeric vector) ~ factor, data = data.frame 형식으로 사용하는 방법이며, 두번째는 x, y numeric vector 를 이용하는 방법입니다. 아래에 2가지 방법 모두를 순서대로 소개하였으며, 분석 결과는 똑같습니다.

(3-1) 두 모집단의 모평균 차이에 대한 추정과 검정 : t-test

- 방법 1 : x ~ factor, data = data.frame

> 
> # independent two sample t-test : t-test()
> # H0 : price of USA vs. non-USA is equal
> # H1 : price of USA vs. non-USA is not equal
> 
> ##-- way 1. x ~ factor, data = data.frame
> t.test(Price ~ Origin, 
+        data = Cars93, 
+        alternative = c("two.sided"), # c("two.sided", "less", "greater")
+        var.equal = FALSE, 
+        conf.level = 0.95)

	Welch Two Sample t-test

data:  Price by Origin
t = -0.95449, df = 77.667, p-value = 0.3428
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
 -5.974255  2.102311
sample estimates:
    mean in group USA mean in group non-USA 
             18.57292              20.50889

(3-2) 두 모집단의 모평균 차이에 대한 추정과 검정 : t-test

- 방법 2 : x, y numeric vectors

> 
> ##--- way 2. x, y numeric vectors

>
> # x, y numeric vector indexing
> Price_USA <- Cars93[which(Cars93$Origin == c("USA")), c("Price")]
> Price_nonUSA <- Cars93[which(Cars93$Origin == c("non-USA")), c("Price")]
>
> # t-test() 
> t.test(Price_USA, Price_nonUSA, 
+        alternative = c("two.sided"),  # c("two.sided", "less", "greater")
+        var.equal = FALSE, 
+        conf.level = 0.95)

	Welch Two Sample t-test

data:  Price_USA and Price_nonUSA
t = -0.95449, df = 77.667, p-value = 0.3428
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
 -5.974255  2.102311
sample estimates:
mean of x mean of y 
 18.57292  20.50889

P-value 가 0.3428로서 5% 유의수준(significance level) 하에서 양측검정(two-sided test) 결과 'USA와 non-USA의 가격의 모평균에는 차이가 없다'는 귀무가설(H0)을 채택(accept)할 수 있게 되었습니다.

그리고, 두 모집단의 모평균 차이에 대한 95% 신뢰구간은 (-5.974255, 2.102311) 로 나왔습니다.

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

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

728x90

저작자표시 비영리 변경금지

'R 분석과 프로그래밍 > R 통계분석' 카테고리의 다른 글

R (3) 두 모집단의 모비율 차이에 대한 추정과 검정 : prop.test() (9)	2015.10.18
R (2) 짝을 이룬 표본에 대한 평균 차이에 대한 추정과 검정(paired sample t-test) : t.test(paired=TRUE) (0)	2015.10.15
R 두 연속형 변수 상관관계 : 공분산(covariance) cov(), 상관계수(correlation coefficient) cor() (4)	2015.10.12
R 연속형 변수 (2) 그룹별(요인별) 요약통계 비교하기 : tapply(var, factor, summary), by(), aggregate(), summaryBy(), describeBy() (4)	2015.10.12
R 연속형 변수 (1) 요약통계 한번에 보기 : summary(), stat.desc(), describe() (2)	2015.10.12

Posted by Rfriend

R 두 연속형 변수 상관관계 : 공분산(covariance) cov(), 상관계수(correlation coefficient) cor()

R 분석과 프로그래밍/R 통계분석 2015. 10. 12. 23:29

다변량 통계분석의 경우 두 변수 간의 상관관계에 아주 많이 의존합니다. 따라서 상관관계 두 연속형 변수 (continuous variable)의 상관관계를 분석하는 것이 매우 중요한데요, 상관관계 분석 기법으로는

통계량 분석

(1) 공분산 (covariance)

(2) 상관계수 (correlation coefficient)

그래프 분석

(1) 산점도 (scatter plot)

(2) 산점도 행렬(scatter matrix plot)

(3) 상관계수행렬(correlation coefficient plot)

등이 있습니다.

그래프는 오른쪽 메뉴의 'R 그래프/시각화'에서 ggplot 을 이용해서 그래프 작성하는 방법을 이미 소개하였으므로, 이번 포스팅에서는 통계량 분석에 대해서만 알아보도록 하겠습니다.

R 분석 예제 데이터는 MASS 패키지 내 Cars93 데이터프레임의 고속도로연비(MPG.highway)와 무게(Weight) 연속형 변수를 사용하겠습니다.

> library(MASS)
> str(Cars93)
'data.frame':	93 obs. of  27 variables:
 $ Manufacturer      : Factor w/ 32 levels "Acura","Audi",..: 1 1 2 2 3 4 4 4 4 5 ...
 $ Model             : Factor w/ 93 levels "100","190E","240",..: 49 56 9 1 6 24 54 74 73 35 ...
 $ Type              : Factor w/ 6 levels "Compact","Large",..: 4 3 1 3 3 3 2 2 3 2 ...
 $ Min.Price         : num  12.9 29.2 25.9 30.8 23.7 14.2 19.9 22.6 26.3 33 ...
 $ Price             : num  15.9 33.9 29.1 37.7 30 15.7 20.8 23.7 26.3 34.7 ...
 $ Max.Price         : num  18.8 38.7 32.3 44.6 36.2 17.3 21.7 24.9 26.3 36.3 ...
 $ MPG.city          : int  25 18 20 19 22 22 19 16 19 16 ...
 $ MPG.highway       : int  31 25 26 26 30 31 28 25 27 25 ...
 $ AirBags           : Factor w/ 3 levels "Driver & Passenger",..: 3 1 2 1 2 2 2 2 2 2 ...
 $ DriveTrain        : Factor w/ 3 levels "4WD","Front",..: 2 2 2 2 3 2 2 3 2 2 ...
 $ Cylinders         : Factor w/ 6 levels "3","4","5","6",..: 2 4 4 4 2 2 4 4 4 5 ...
 $ EngineSize        : num  1.8 3.2 2.8 2.8 3.5 2.2 3.8 5.7 3.8 4.9 ...
 $ Horsepower        : int  140 200 172 172 208 110 170 180 170 200 ...
 $ RPM               : int  6300 5500 5500 5500 5700 5200 4800 4000 4800 4100 ...
 $ Rev.per.mile      : int  2890 2335 2280 2535 2545 2565 1570 1320 1690 1510 ...
 $ Man.trans.avail   : Factor w/ 2 levels "No","Yes": 2 2 2 2 2 1 1 1 1 1 ...
 $ Fuel.tank.capacity: num  13.2 18 16.9 21.1 21.1 16.4 18 23 18.8 18 ...
 $ Passengers        : int  5 5 5 6 4 6 6 6 5 6 ...
 $ Length            : int  177 195 180 193 186 189 200 216 198 206 ...
 $ Wheelbase         : int  102 115 102 106 109 105 111 116 108 114 ...
 $ Width             : int  68 71 67 70 69 69 74 78 73 73 ...
 $ Turn.circle       : int  37 38 37 37 39 41 42 45 41 43 ...
 $ Rear.seat.room    : num  26.5 30 28 31 27 28 30.5 30.5 26.5 35 ...
 $ Luggage.room      : int  11 15 14 17 13 16 17 21 14 18 ...
 $ Weight            : int  2705 3560 3375 3405 3640 2880 3470 4105 3495 3620 ...
 $ Origin            : Factor w/ 2 levels "USA","non-USA": 2 2 2 2 2 1 1 1 1 1 ...
 $ Make              : Factor w/ 93 levels "Acura Integra",..: 1 2 4 3 5 6 7 9 8 10 ...

상식적으로 생각해보면 차 무게가 많이 나갈수록 고속도로연비가 안좋게 나올텐데요, 아래의 산점도를 보면 역(-)의 상관관계가 있음을 알 수 있습니다.

> # Scatter Plot of MPG.highway ~ Weight with linear regression line

> plot(MPG.highway ~ Weight, data = Cars93, 
+      main = "Scatter Plot of MPG.highway ~ Weight", 
+      xlab = "Weight", ylab = "MPG.highway")
> abline(lm(MPG.highway ~ Weight, data = Cars93), col = "blue", lty = 2)

그럼, 먼저 공분산(covariance)을 알아보고, 그 다음으로 상관계수(correlation coefficient)를 알아보겠습니다.

(1) 공분산 (covariance) : cov()

공분산의 이해를 돕기 위해 x가 증가할 때 y도 역시 증가하는 경우를 가정해보면서 설명을 해보겠습니다.

- μ(x) 보다 큰 x값에 대응하는 y값도 또한 μ(y)보다 크므로 (X-μ(x)(Y-μ(y))는 양수가 됩니다.(+*+=+)

- μ(x) 보다 작은 x값에 대응하는 y값도 또한 μ(y)보다 작으므로 (X-μ(x)(Y-μ(y))는 양수가 됩니다.(-*-=+)

반대로, x가 증가할 때 y는 감소하는 경우를 가정해봅시다. 이 경우에는,

- μ(x) 보다 큰 x값에 대응하는 y값은 μ(y)보다 작으므로 (X-μ(x)(Y-μ(y))는 음수가 됩니다.(+*-=-)

- μ(x) 보다 작은 x값에 대응하는 y값은 μ(y)보다 크므로 (X-μ(x)(Y-μ(y))는 음수가 됩니다.(-*+=-)

즉, 공분산(covriance) 값이 '+'이면 x와 y 두 변수는 같은 방향으로 움직이며 (x가 증가하면 y도 증가, x가 감소하면 y도 감소), 공분산이 '-'이면 x와 y 두 변수는 서로 다른 방향으로 움직인다는 것을 알 수 있습니다.

R의 cov() 함수를 사용해서 고속도로연비(MPG.highway)와 무게(Weight)와의 공분산을 2가지 방법을 사용해서 구해보겠습니다.

> ##-- covariance 
> # way 1. cov(x, y)
> with(Cars93, cov(x=MPG.highway, 
+                  y=Weight, 
+                  use="complete.obs", 
+                  method=c("pearson")))
[1] -2549.655
> # way 2. cov(x, y=NULL)
> x <- Cars93[c("MPG.highway", "Weight")]
> with(Cars93, cov(x, 
+                  y=NULL, 
+                  use="complete.obs", 
+                  method=c("pearson")))
            MPG.highway     Weight
MPG.highway     28.4273  -2549.655
Weight       -2549.6546 347977.893

way 1. cov(x, y)를 사용(x, y는 벡터)하든지, way 2. cov(x, y=NULL)을 사용(x는 데이터프레임)하든지 결과는 pearson covariance 가 -2549.6546 으로 음의 상관관계를 가지는 것으로 나왔습니다.

그런데 단위가 -2549.6546 이면 어느정도로 음의 상관관계라고 하는 것인지 감이 안잡힙니다. 이때 공분산을 표준화한 상관계수를 사용하면 두 변수간의 상관관계의 정도를 가늠할 수 있게 됩니다.

(2) 상관계수 (pearson correlation coefficient) : cor()

두 변수간의 관련성을 나타내는 측도서 상관계수는 다음의 특징을 가집니다.

1) 표준화변수들의 공분산은 상관계수가 된다.

2) 상관계수는 -1 ≤ ρ ≤ 1 의 범위의 값을 가진다.

3) 상관계수 ρ가 +1에 가까울수록 강한 양의 선형관계를 가지며,

-1에 가까울수록 강한 음의 선형관계를 가진다.

4) 상관계수 ρ = 0 이면 두 변수간의 선형관계는 없으며,

0에 가까울수록 선형관계가 약해진다.

(단, 비선형 관계를 가질 수는 있음. 그래프 분석 병행 필요)

5) 위치 변환이나 척도 변환 후에도 상관계수는 변함이 없다.

R의 cor() 함수를 사용해서 고속도로연비(MPG.highway)와 무게(Weight)의 pearson 상관계수(pearson correlation coefficient)를 구해보면 다음과 같습니다.

> ##-- correlation coefficient
> # way 1. cor(x, y)
> with(Cars93, cor(x=MPG.highway, 
+                  y=Weight, 
+                  use="complete.obs", 
+                  method=c("pearson")))
[1] -0.8106581
> 
> 
> x <- Cars93[c("MPG.highway", "Weight")]
> cor(x, 
+     y=NULL, 
+     use="complete.obs", 
+     method=c("pearson"))
            MPG.highway     Weight
MPG.highway   1.0000000 -0.8106581
Weight       -0.8106581  1.0000000
> 
> 
> # options
> # use = c("everythig", "all.obs", "complete.obs", "na.or.complete", pairwise.complete.obs")
> # method = c("pearson", "kendall", "spearman")

고속도로연비(MPG.highway)와 무게(Weight)의 pearson 상관계수는 -0.8106581 으로서 매우 강한 음의 상관관계를 가지는 것으로 나왔습니다.

cov()와 cor() 함수에 공통으로 use와 method option이 사용되었는데요, 이에 대해 추가 설명을 하겠습니다.

구분	option	내용
use = (결측값 처리)	"everything" (default)	결측값이 있을 경우 NA로 계산 결과 제시
	"all.obs"	결측값이 있을 경우 오류 발생 ("Error in cor(x = a, y = b, use = "all.obs", method = c("pearson")) : cov/cor에 결측치들이 있습니다 " error message)
	"complete.obs"	결측값이 있는 case는 모두 제거된 상태에서 상관계수 계산
	"pairwise.complete.obs"	상관계수가 계산되는 변수들만을 대상으로 결측값이 있는 case 제거한 상관계수 계산
method = (상관계수 통계량)	"pearson" (default)	Pearson correlation coefficient 지정, 가장 일반적으로 사용
	"kendall"	Kendall의 순위상관계수 혹은 Kendall의 τ (tau) 지정, 비모수 상관계수 계산 (정규성 불충족 시)
	"spearman"	Spearman의 순위상관계수 혹은 Spearman이 ρ(rho) 지정, 비모수 상관계수 계산 (정규성 불충족 시)

* kendall τ , Spearman ρ 수식은 생략함

아래에는 결측값(NA, missing value)이 포함되어 있을 경우 공분산이나 상관계수 계산 시에 결측값을 처리하는 옵션인 'use' 의 세부 옵션별("everything", "all.obs", "complete.obs", "pairwise.complete.obs")로 어떻게 다른지를 간단한 예제를 들어보았습니다.

a 벡터의 5번째 원소에 NA가 들어있고, c 벡터의 1번째 원소에 NA가 들어있는 예제입니다. use 옵션을 좀더 명확히 이해하는데 도움이 되기를 바랍니다.

- use = "everything" : NA가 포함되어 있으므로 결과는 "NA"

- use = "all.obs" : NA가 포함되어 있으므로 다음과 같은 error message => "cov/cor에 결측치들이 있습니다"

- use = "complete.obs" : a, b, c 벡터에서 한개라도 NA가 있으면 그 순번째의 a, b, c 벡터의 원소(만약 데이터 프레임이라면 해당 row 전체, 즉 해당 observation case의 모든 변수)는 통째로 삭제하라는 뜻이므로 a, b, c 벡터 모두의 1번째와 5번째 원소 삭제한 후에 상관계수 계산

- use = "pairwise.compelte.obs" : a와 b의 짝 (pairwise b/w a and b)에서는 5번째 원소 삭제, a와 c 의 짝(pairwise b/w a and d)에서는 1번째와 5번째 원소 삭제, b와 c의 짝(pairwise b/w b and c)에서는 1번째 원소 삭제한 후에 상관계수 계산

> ##------------
> ## cov(), cor()
> ## use : an optional character string giving a method for computing covariances 
> ##       in the presence of missing values
> ##------------
> 
> # dataset for correlation analysis with NAs
> a <- c(10, 9, 8, 6, NA)
> b <- c(9, 6, 4, 3, 2)
> c <- c(NA, 13, 18, 20, 25)
> 
> x <- data.frame(a, b, c)
> 
> 
> ##------------
> # use = "everything"  
> #  ; (default setting) NAs will propagate conceptually, 
> #    i.e., a resulting value will be NA whenever one of its contributing observations is NA
> cor(x, 
+     y = NULL, 
+     use = "everything", 
+     method = c("pearson"))
   a  b  c
a  1 NA NA
b NA  1 NA
c NA NA  1
> 
> 
> ##------------
> # use = "all.obs" 
> # ; the presence of missing observations will produce an error
> cor(x, 
+     y = NULL, 
+     use = "all.obs", 
+     method = c("pearson"))
Error in cor(x, y = NULL, use = "all.obs", method = c("pearson")) : 
  cov/cor에 결측치들이 있습니다
> 
> 
> ##------------
> # use = "complete.obs" 
> #  ; missing values are handled by casewise deletion 
> #    (and if there are no complete cases, that gives an error)
> cor(x, 
+     y = NULL, 
+     use = "complete.obs", 
+     method = c("pearson"))
           a          b          c
a  1.0000000  0.9285714 -0.9078413
b  0.9285714  1.0000000 -0.9986254
c -0.9078413 -0.9986254  1.0000000
> 
> 
> # double check
> a <- c(10, 9, 8, 6, NA)
> b <- c(9, 6, 4, 3, 2)
> c <- c(NA, 13, 18, 20, 25)
> 
>   # use = "complete.obs" : delete 1st and 5th elements of a, b, c vectors
> a_complete.obs <- c(9, 8, 6)
> b_complete.obs <- c(6, 4, 3)
> c_complete.obs <- c(13, 18, 20)
> x_complete.obs <- data.frame(a_complete.obs, b_complete.obs, c_complete.obs)
> 
> cor(x_complete.obs, 
+     y = NULL, 
+     use = "all.obs", # because it is calculated after deleting 1st and 5th elements
+     method = c("pearson"))
               a_complete.obs b_complete.obs c_complete.obs
a_complete.obs      1.0000000      0.9285714     -0.9078413
b_complete.obs      0.9285714      1.0000000     -0.9986254
c_complete.obs     -0.9078413     -0.9986254      1.0000000
> 
> 
> ##------------
> # use = "pairwise.complete.obs" 
> #  ; the correlation or covariance between each pair of 
> #    variables is computed using all complete pairs of observations on those variables
> #    only works with the "pearson" method
> cor(x, 
+     y = NULL, 
+     use = "pairwise.complete.obs", 
+     method = c("pearson"))
           a          b          c
a  1.0000000  0.9221389 -0.9078413
b  0.9221389  1.0000000 -0.9824719
c -0.9078413 -0.9824719  1.0000000
> 
> # double check 
> a <- c(10, 9, 8, 6, NA)
> b <- c(9, 6, 4, 3, 2)
> c <- c(NA, 13, 18, 20, 25)
>   
>   # use = "pairwise.complete.obs" : delete 5th element of a and b pair vectors
> a_ab_pairwise.complete.obs <- c(10, 9, 8, 6)
> b_ab_pairwise.complete.obs <- c(9, 6, 4, 3)
> 
> cor(x = a_ab_pairwise.complete.obs, 
+     y = b_ab_pairwise.complete.obs, 
+     use = "all.obs", # because it is calculated after deleting 5th element
+     method = c("pearson"))
[1] 0.9221389
> 
>   # use = "pairwise.complete.obs" : delete 1st and 5th elements of a and c pair vectors
> a_ac_pairwise.complete.obs <- c(9, 8, 6)
> c_ac_pairwise.complete.obs <- c(13, 18, 20)
> 
> cor(x = a_ac_pairwise.complete.obs, 
+     y = c_ac_pairwise.complete.obs, 
+     use = "all.obs", 
+     method = c("pearson"))
[1] -0.9078413
> 
>   # use = "pairwise.complete.obs" : delete 1st element of b and c pair vectors
> b_bc_pairwise.complete.obs <- c(6, 4, 3, 2)
> c_bc_pairwise.complete.obs <- c(13, 18, 20, 25)
> 
> cor(x = b_bc_pairwise.complete.obs, 
+     y = c_bc_pairwise.complete.obs, 
+     use = "all.obs", 
+     method = c("pearson"))
[1] -0.9824719

(3) 여러개의 숫자형 변수 간 상관계수 결과를 세로로 긴 형태(Long Format)로 저장하기

여러개의 숫자형 변수로 이루어진 데이터 프레임에서 숫자형 변수의 짝(pairs) 별로 상관계수를 구해서 세로로 길쭉한 형태(long format)으로 결과를 저장하고 싶다면 아래의 for loop 문을 사용한 예제를 참고하세요.

표준정규분포로 부터 10개 칼럼, 5개 관측치의 난수를 생성하여 예제 숫자형 데이터 프레임을 만들었습니다.

for loop 순환문으로 칼럼 갯수만큼 i, j 를 돌리면 두개의 변수 쌍(pairs)이 두 번씩 출현할텐데요, 두 변수의 순서가 바뀌어도 상관계수 값은 동일하므로 두 변수 쌍은 한번만 계산하도록 for loop 문을 작성했습니다. 그리고 자기 자신과의 상관계수는 항상 '1'이므로 제외하도록 했습니다.

# generate 50 random numbers from standard normal distribution X~N(0, 1)

set.seed(1004)

x <- rnorm(50, mean=0, sd=1)

df <- data.frame(matrix(x, nrow=5, ncol=10, byrow=TRUE))

# X1 X2 X3 X4 X5 X6 X7 X8 X9 X10

# 1 -0.6078408 0.7676292 -0.1642718 -0.02877676 0.01355787 -0.7157169 -1.477976 0.8280146 -0.8951765 1.8841561

# 2 -1.8952567 1.1259354 0.1494180 -0.21511242 2.63756778 -0.2678609 1.123417 0.4842308 0.6318069 -0.8180967

# 3 -0.7470136 -0.5132835 0.8709808 0.83062565 -1.03451768 -1.6140827 2.169269 0.6719930 0.4054343 -0.8316793

# 4 0.9159031 0.3905532 -0.6347263 -1.00658131 0.35622463 0.1083284 -1.712486 0.6846194 -0.6083049 -0.6192043

# 5 -0.7671233 0.3633518 -0.7325597 1.23058116 0.09165232 1.5699387 1.309473 0.9069780 -1.3729150 -0.1163513

# pearson's correlation coefficients in a long format

cor_long <- data.frame() # blank dataframe to save the corr results later

col_nm <- colnames(df)

n <- 1

for (i in 1:(ncol(df)-1)){

for (j in (i+1):ncol(df)) {

cor_long[n, 1] <- col_nm[i]

cor_long[n, 2] <- col_nm[j]

cor_long[n, 3] <- cor(df[,i], df[,j])

n <- n+1

}

# rename column names

colnames(cor_long) <- c("col_1", "col_2", "corr_coef")

cor_long

# col_1 col_2 corr_coef

# 1 X1 X2 -0.332162830

# 2 X1 X3 -0.448170165

# 3 X1 X4 -0.438930937

# 4 X1 X5 -0.497470094

# 5 X1 X6 0.100596534

# 6 X1 X7 -0.643972517

# 7 X1 X8 0.346978141

# 8 X1 X9 -0.465219120

# 9 X1 X10 0.051698753

# 10 X2 X3 -0.440364182

# 11 X2 X4 -0.454194089

# 12 X2 X5 0.847973001

# 13 X2 X6 0.342903949

# 14 X2 X7 -0.413224712

# 15 X2 X8 -0.235257160

# 16 X2 X9 -0.064519846

# 17 X2 X10 0.314659871

# 18 X3 X4 0.224054384

# 19 X3 X5 -0.167990551

# 20 X3 X6 -0.865225316

# 21 X3 X7 0.560975754

# 22 X3 X8 -0.500395969

# 23 X3 X9 0.805120716

# 24 X3 X10 -0.244220838

# 25 X4 X5 -0.419495665

# 26 X4 X6 0.170965784

# 27 X4 X7 0.755460153

# 28 X4 X8 0.482063832

# 29 X4 X9 -0.181887756

# 30 X4 X10 -0.001894032

# 31 X5 X6 0.239111079

# 32 X5 X7 -0.039280725

# 33 X5 X8 -0.632907538

# 34 X5 X9 0.356088708

# 35 X5 X10 -0.209073460

# 36 X6 X7 -0.074444629

# 37 X6 X8 0.457440745

# 38 X6 X9 -0.681812940

# 39 X6 X10 -0.023974644

# 40 X7 X8 -0.189050606

# 41 X7 X9 0.442853344

# 42 X7 X10 -0.552145749

# 43 X8 X9 -0.913090461

# 44 X8 X10 0.582373314

# 45 X9 X10 -0.556365998

R로 상관계수 행렬 그림을 그리는 방법은 https://rfriend.tistory.com/82 를 참고하세요.

두 변수 x, y별 상관계수를 산점도로 나타낸 아래의 그림을 참조하시기 바랍니다. 상관계수가 '0'이면 '선형관계'는 없는 것이지만 '비선형관계'는 있을 수도 있음에 유의하시기 바랍니다. 두번째 그림에서 보면 상관계수가 '1'인 것은 x, y가 선으로 되어있는 관계이며, 기울기는 상관없음을 알 수 있습니다. y 값이 '0'인 선은 상관계수 값이 없습니다.

[ Examples of correlation coefficients ]

* source : https://en.wikipedia.org/wiki/Pearson_product-moment_correlation_coefficient

상관관계(correlation) 분석할 때 인과관계(causality)와 혼동하는 경우를 종종 봅니다. 상관관계(correlation)는 단순히 두 변수 x와 y간의 상관성만을 보는 것이며, x와 y가 상관계수가 높다고 해서 x가 y의 원인이라거나 y가 x의 원인이라고 해석하면 안됩니다.

인과관계(causality)에는 아래처럼 6가지 종류의 인과관계가 있으며, 시간의 흐름(x chronologically precedes y), 영향을 끼치는 방향(x causally precedes y)을 함게 고려하여야 하므로 사용, 해석에 주의를 요합니다.

[ 6 different kinds of causality ]

There are six different kinds of causality within a model.

1. A direct causal relationship is one in which a variable, X, is a direct cause of another variable, Y (i.e. it is the immediate determinant of Y within the context of the theoretical system).

2. An indirect causal relationship is one in which X exerts a causal impact on Y, but only through its impact on a third variable, Z.

3. A spurious relationship is one in which X and Y are related, but only because of a common cause, Z. There is no formal causal link between X and Y.

4. A bi-directional or reciprocal causal relationship is one in which X has a causal influence on Y, which in turn, has a causal impact on X.

5. An unanalyzed relationship is one in which X and Y or related, but the source of the relationship is unspecified.

6. A moderated causal relationship is one in which the relationship between X ad Y is moderated by a third variable. In other words, the nature of the relationship between X and Y varies, depending on the value of Z.

From Jaccard, J., & Turrisi, R. (2003). Interaction Effects in Multiple Regression. Thousand Oaks, CA: Sage.

* source : http://depts.washington.edu/methods/causality.html

위의 영어 설명을 예를 들어서 도식화하면 아래와 같습니다.

[ 인과관계의 종류 및 예시 ]

산점도, 산점도 행렬, 상관계수 행렬 그림 그리기는 아래 포스팅을 참조하시기 바랍니다.

☞ ggplot2의 geom_point() 함수를 이용한 한개 변수의 산점도 그리기

☞ R Graphics 패키지 pairs() 함수를 사용한 산점도 행렬 그리기

☞ corrplot 패키지를 이용한 상관계수 행렬 그림 그리기

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

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

728x90

저작자표시 비영리 변경금지

'R 분석과 프로그래밍 > R 통계분석' 카테고리의 다른 글

R (2) 짝을 이룬 표본에 대한 평균 차이에 대한 추정과 검정(paired sample t-test) : t.test(paired=TRUE) (0)	2015.10.15
R (1) 독립된 두 표본의 모평균 차이에 대한 추정과 검정 : t.test() (11)	2015.10.15
R 연속형 변수 (2) 그룹별(요인별) 요약통계 비교하기 : tapply(var, factor, summary), by(), aggregate(), summaryBy(), describeBy() (4)	2015.10.12
R 연속형 변수 (1) 요약통계 한번에 보기 : summary(), stat.desc(), describe() (2)	2015.10.12
R 일변량 연속형 자료 요약 - (3) 분포형태와 대칭정도(distrubution, symmetry) : 왜도(skewness), 첨도(kurtosis), 분위수-분위수(Quantile-Qunantile) (0)	2015.10.12

Posted by Rfriend

R 연속형 변수 (2) 그룹별(요인별) 요약통계 비교하기 : tapply(var, factor, summary), by(), aggregate(), summaryBy(), describeBy()

R 분석과 프로그래밍/R 통계분석 2015. 10. 12. 18:20

이전 포스팅에서 연속형 변수(continuous variable)의 중심화 경향(central tendancy), 퍼짐 정도(dispersion), 분포형태 및 대칭정도 (distribution and symmetry) 등의 R 함수에 대해 알아보았습니다.

R에서는 이들 개별 통계량에 대해 개별 R 함수를 제공함과 동시에, 참 편리하게도 한방에(!) 연속형 변수의 기술통계량을 볼 수 있게 해주는 함수들이 있는데요,

(1) 연속형 변수 요약통계 한번에 보기

: summary(), stat.desc(), describe()

(2) 연속형 변수 그룹별(요인별) 요약통계 비교하기

: tapply(var, factor, summary), aggregate(), summaryBy(), describe.by()

중에서 이번에는 R을 이용한 (2) 연속형 변수 그룹별(요인별) 요약통계 비교하기를 소개하겠습니다.

[ 연속형 변수 요약통계 한번에 보기 package & function 요약 ]

category	package	function	statistics
(1) descriptive statistics summary	base	summary()	min, 1Q, median, mean, 3Q, max
	pastecs	stat.desc()	nbr.val, nbr.null, nbr.na, min, max, range, sum, median, mean, var, std.dev, coef.var, skewness, kurtosis, normtest.W, normtest.p, SE.mean, CI.mean.p
	psych	describe()	n, mean, std.dev, median, trimmed, mad, min, max, range, skew, kurtosis, se
(2) descriptive statistics comparison by group (factor)	base	tapply(var, factor, summary)	summary: min, 1Q, median, mean, 3Q, max, (user defined) functions...
	base	by()	(user defined) functions...
	stats	aggregate()	(user defined) functions...
	doBy	summaryBy()	(user defined) functions...
	psych	describeBy()	n, mean, sd, median, trimmed, mad, min, max, range, skew, kurtosis, se

MASS package 내 Cars93 데이터프레임의 차종(Type)별 가격(Price), 고속도로연비(MPG.highway) 연속형 변수에 대해서 요약통계를 계산해보겠습니다.

> library(MASS)
> str(Cars93)
'data.frame':	93 obs. of  27 variables:
 $ Manufacturer      : Factor w/ 32 levels "Acura","Audi",..: 1 1 2 2 3 4 4 4 4 5 ...
 $ Model             : Factor w/ 93 levels "100","190E","240",..: 49 56 9 1 6 24 54 74 73 35 ...
 $ Type              : Factor w/ 6 levels "Compact","Large",..: 4 3 1 3 3 3 2 2 3 2 ...
 $ Min.Price         : num  12.9 29.2 25.9 30.8 23.7 14.2 19.9 22.6 26.3 33 ...
 $ Price             : num  15.9 33.9 29.1 37.7 30 15.7 20.8 23.7 26.3 34.7 ...
 $ Max.Price         : num  18.8 38.7 32.3 44.6 36.2 17.3 21.7 24.9 26.3 36.3 ...
 $ MPG.city          : int  25 18 20 19 22 22 19 16 19 16 ...
 $ MPG.highway       : int  31 25 26 26 30 31 28 25 27 25 ...
 $ AirBags           : Factor w/ 3 levels "Driver & Passenger",..: 3 1 2 1 2 2 2 2 2 2 ...
 $ DriveTrain        : Factor w/ 3 levels "4WD","Front",..: 2 2 2 2 3 2 2 3 2 2 ...
 $ Cylinders         : Factor w/ 6 levels "3","4","5","6",..: 2 4 4 4 2 2 4 4 4 5 ...
 $ EngineSize        : num  1.8 3.2 2.8 2.8 3.5 2.2 3.8 5.7 3.8 4.9 ...
 $ Horsepower        : int  140 200 172 172 208 110 170 180 170 200 ...
 $ RPM               : int  6300 5500 5500 5500 5700 5200 4800 4000 4800 4100 ...
 $ Rev.per.mile      : int  2890 2335 2280 2535 2545 2565 1570 1320 1690 1510 ...
 $ Man.trans.avail   : Factor w/ 2 levels "No","Yes": 2 2 2 2 2 1 1 1 1 1 ...
 $ Fuel.tank.capacity: num  13.2 18 16.9 21.1 21.1 16.4 18 23 18.8 18 ...
 $ Passengers        : int  5 5 5 6 4 6 6 6 5 6 ...
 $ Length            : int  177 195 180 193 186 189 200 216 198 206 ...
 $ Wheelbase         : int  102 115 102 106 109 105 111 116 108 114 ...
 $ Width             : int  68 71 67 70 69 69 74 78 73 73 ...
 $ Turn.circle       : int  37 38 37 37 39 41 42 45 41 43 ...
 $ Rear.seat.room    : num  26.5 30 28 31 27 28 30.5 30.5 26.5 35 ...
 $ Luggage.room      : int  11 15 14 17 13 16 17 21 14 18 ...
 $ Weight            : int  2705 3560 3375 3405 3640 2880 3470 4105 3495 3620 ...
 $ Origin            : Factor w/ 2 levels "USA","non-USA": 2 2 2 2 2 1 1 1 1 1 ...
 $ Make              : Factor w/ 93 levels "Acura Integra",..: 1 2 4 3 5 6 7 9 8 10 ...

각 package별 R 함수 사용예는 아래와 같습니다.

(1) base package : tapply(var, factor, summary)

> ##-- (1) base package : with(dataset, tapply(variable, factor, summary))
> 
> with(Cars93, tapply(Price, Type, summary))
$Compact
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
     11      13      16      18      21      32 

$Large
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
     18      20      21      24      27      36 

$Midsize
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
     14      17      26      27      34      62 

$Small
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
    7.4     8.6    10.0    10.2    11.3    15.9 

$Sporty
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
     10      14      17      19      22      38 

$Van
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
     16      19      19      19      20      23

> with(Cars93, tapply(MPG.highway, Type, summary))
$Compact
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
     26      28      30      30      31      36 

$Large
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
     25      26      26      27      28      28 

$Midsize
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
     22      25      26      27      29      31 

$Small
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
     29      33      33      36      37      50 

$Sporty
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
     24      25      28      29      31      36 

$Van
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
     20      21      22      22      23      24

tapply함수를 사용해서 특정 변수에 요인(factor)별 summary()함수를 적용하였습니다. 가격(Price)과 고속도로연비(MPG.highway) 를 한꺼번에 못하고 두번에 나누어서 적용해야 하며, 제시된 포맷도 일괄 비교해서 보기에는 조금 불편한 점이 있습니다.

(2) base package : by()

> ##-- (2) base package : by()
> fun_summary <- function(x, ...) {
+   c(n=sum(!is.na(x)), mean=mean(x, ...), sd=sd(x, ...))
+ }
> 
> by(Cars93[c("Price", "MPG.highway")], # dataset$variable
+    Cars93$Type, # by Group(Factor)
+    function(x) sapply(x, fun_summary, na.rm=T)) # function
Cars93$Type: Compact
     Price MPG.highway
n     16.0        16.0
mean  18.2        29.9
sd     6.7         2.9
-------------------------------------------------------------------------- 
Cars93$Type: Large
     Price MPG.highway
n     11.0        11.0
mean  24.3        26.7
sd     6.3         1.3
-------------------------------------------------------------------------- 
Cars93$Type: Midsize
     Price MPG.highway
n       22        22.0
mean    27        26.7
sd      12         2.5
-------------------------------------------------------------------------- 
Cars93$Type: Small
     Price MPG.highway
n       21        21.0
mean    10        35.5
sd       2         5.6
-------------------------------------------------------------------------- 
Cars93$Type: Sporty
     Price MPG.highway
n       14        14.0
mean    19        28.8
sd       8         3.6
-------------------------------------------------------------------------- 
Cars93$Type: Van
     Price MPG.highway
n      9.0         9.0
mean  19.1        21.9
sd     1.9         1.5

두번째는 base package의 by()함수를 사용하는 방법인데요, 위의 예에서는 관측값 개수 n, 평균 mean, 표준편차 sd 를 그룹(요인)별로 일괄 계산하는 script가 되겠습니다.

by() 함수는 두 개 이상의 기술통계량을 그룹(요인)별로 계산하고 싶다면 사용자 정의 함수(user defined function)를 미리 만들어놓아야 하며, by() 함수 내에 sapply()를 적용해야 하므로 초보자가 이용하기에는 어려운 방법이 되겠습니다. 너무 복잡해서 외우기도 사실상 어렵습니다. -_-;

사용자 정의함수에서 '...' 은 '결측값 제거 후 계산'과 같은 추가 옵션을 설정할 수 있도록 하는 명령문입니다.

(3) stats package : aggregate()

> ##-- (3) stats package : aggregate()
> library(stats)
> fun_summary_2 <- function(x, ...) {
+   c(n=sum(!is.na(x)), mean=mean(x, ...), sd=sd(x, ...))
+ }
> 
> aggregate(Cars93[c("Price", "MPG.highway")], # dataset$variable
+           by=list(Car_Type=Cars93$Type), # by Group(Factor)
+           fun_summary_2, na.rm=T) # function
  Car_Type Price.n Price.mean Price.sd MPG.highway.n MPG.highway.mean MPG.highway.sd
1  Compact    16.0       18.2      6.7          16.0             29.9            2.9
2    Large    11.0       24.3      6.3          11.0             26.7            1.3
3  Midsize    22.0       27.2     12.3          22.0             26.7            2.5
4    Small    21.0       10.2      2.0          21.0             35.5            5.6
5   Sporty    14.0       19.4      8.0          14.0             28.8            3.6
6      Van     9.0       19.1      1.9           9.0             21.9            1.5

세번째로, stats package의 aggregate() 함수는 그룹(요인) 자리에 list() 로 들어가며, by()함수처럼 두개 이상의 통계량을 보려면 사용자 정의 함수를 미리 만들어서 사용해야 하므로 초보자에게 역시 쉽지 않습니다.

통계량들이 제시되는 포맷이 tapply(summary), by()와는 달리, 위에서 보는 것처럼 row에 그룹(요인)이 위치하고 column에 각 통계량들이 나열되어 있습니다.

(4) doBy package : summaryBy()

> ##-- (4) doBy package : summaryBy()
> install.packages("doBy")
Installing package into ‘C:/Users/user/Documents/R/win-library/3.2’
(as ‘lib’ is unspecified)
trying URL 'http://cran.rstudio.com/bin/windows/contrib/3.2/doBy_4.5-13.zip'
Content type 'application/zip' length 3431577 bytes (3.3 MB)
downloaded 3.3 MB

package ‘doBy’ successfully unpacked and MD5 sums checked

The downloaded binary packages are in
	C:\Users\user\AppData\Local\Temp\RtmpqgHvNR\downloaded_packages
> library(doBy)
필요한 패키지를 로딩중입니다: survival
>

> 
> summaryBy(Price + MPG.highway ~ Type, data=Cars93, FUN=c(mean, sd))


     Type Price.mean MPG.highway.mean Price.sd MPG.highway.sd
1 Compact         18               30      6.7            2.9
2   Large         24               27      6.3            1.3
3 Midsize         27               27     12.3            2.5
4   Small         10               35      2.0            5.6
5  Sporty         19               29      8.0            3.6
6     Van         19               22      1.9            1.5

네번째로 doBy package의 summaryBy() 함수는 기본 표현법이 summaryBy( var1 + var2 ~ factor, data, FUN=c(,,)) 으로서 초보자도 사용하기에 매우 쉽고 직관적으로 되어있습니다. 분석결과가 제시되는 포맷도 row에 그룹(요인)이 오고 column에 기초통계량이 와서 일목요연하게 그룹 간 비교하기가 편하게 되어있어서 제가 제일 선호하는 함수입니다.

(5) psych package : describe.by()

> ##-- (5) psych package : describeBy()
> library(psych)
> 
> describeBy(Cars93[c("Price", "MPG.highway")], Cars93$Type, mat=FALSE) # list
group: Compact
            vars  n mean  sd median trimmed mad min max range skew kurtosis   se
Price          1 16   18 6.7     16      18 5.3  11  32    21 0.78    -0.80 1.67
MPG.highway    2 16   30 2.9     30      30 3.0  26  36    10 0.48    -0.82 0.74
-------------------------------------------------------------------------- 
group: Large
            vars  n mean  sd median trimmed mad min max range  skew kurtosis   se
Price          1 11   24 6.3     21      24 3.7  18  36    18  0.82     -1.0 1.91
MPG.highway    2 11   27 1.3     26      27 1.5  25  28     3 -0.07     -1.9 0.38
-------------------------------------------------------------------------- 
group: Midsize
            vars  n mean   sd median trimmed  mad min max range skew kurtosis   se
Price          1 22   27 12.3     26      26 12.7  14  62    48 1.06     0.65 2.61
MPG.highway    2 22   27  2.5     26      27  2.2  22  31     9 0.11    -0.96 0.54
-------------------------------------------------------------------------- 
group: Small
            vars  n mean  sd median trimmed mad  min max range skew kurtosis   se
Price          1 21   10 1.9     10      10 2.1  7.4  16   8.5    1     1.08 0.43
MPG.highway    2 21   35 5.6     33      35 4.4 29.0  50  21.0    1     0.22 1.22
-------------------------------------------------------------------------- 
group: Sporty
            vars  n mean  sd median trimmed mad min max range skew kurtosis   se
Price          1 14   19 8.0     17      19 4.3  10  38    28  1.0     -0.1 2.13
MPG.highway    2 14   29 3.6     28      29 4.4  24  36    12  0.4     -1.1 0.97
-------------------------------------------------------------------------- 
group: Van
            vars n mean  sd median trimmed  mad min max range  skew kurtosis   se
Price          1 9   19 1.9     19      19 0.89  16  23   6.4  0.17     -0.6 0.63
MPG.highway    2 9   22 1.4     22      22 1.48  20  24   4.0 -0.05     -1.7 0.48
> 
> 
> describeBy(Cars93[c("Price", "MPG.highway")], Cars93$Type, mat=TRUE) #matrix
             item  group1 vars  n mean   sd median trimmed   mad  min max range   skew kurtosis
Price1          1 Compact    1 16   18  6.7     16      18  5.34 11.1  32  20.8  0.780    -0.80
Price2          2   Large    1 11   24  6.3     21      24  3.71 18.4  36  17.7  0.815    -1.03
Price3          3 Midsize    1 22   27 12.3     26      26 12.68 13.9  62  48.0  1.060     0.65
Price4          4   Small    1 21   10  2.0     10      10  2.08  7.4  16   8.5  1.007     1.08
Price5          5  Sporty    1 14   19  8.0     17      19  4.30 10.0  38  28.0  1.029    -0.10
Price6          6     Van    1  9   19  1.9     19      19  0.89 16.3  23   6.4  0.166    -0.60
MPG.highway1    7 Compact    2 16   30  2.9     30      30  2.97 26.0  36  10.0  0.483    -0.82
MPG.highway2    8   Large    2 11   27  1.3     26      27  1.48 25.0  28   3.0 -0.068    -1.89
MPG.highway3    9 Midsize    2 22   27  2.5     26      27  2.22 22.0  31   9.0  0.106    -0.96
MPG.highway4   10   Small    2 21   35  5.6     33      35  4.45 29.0  50  21.0  1.021     0.22
MPG.highway5   11  Sporty    2 14   29  3.6     28      29  4.45 24.0  36  12.0  0.399    -1.06
MPG.highway6   12     Van    2  9   22  1.5     22      22  1.48 20.0  24   4.0 -0.049    -1.72
               se
Price1       1.67
Price2       1.91
Price3       2.61
Price4       0.43
Price5       2.13
Price6       0.63
MPG.highway1 0.74
MPG.highway2 0.38
MPG.highway3 0.54
MPG.highway4 1.22
MPG.highway5 0.97
MPG.highway6 0.48

마지막으로 psych package의 describeBy() 함수의 경우 분석가가 임으로 요약통계량을 지정할 수 없게 되어있으며, describe()함수에서 제시되었던 요약통계량들(n, mean, sd, median, trimmed mean, mad, min, max, range, skewness, kurtosis, se 등)이 자동으로 그룹별로 제시됩니다.

mat=FALSE 옵션을 설정하면 list 형식으로 결과가 제시되며, mat=TRUE 옵션을 지정하면 matrix 형식으로 결과가 나타납니다.

위의 package별 function별 제공하는 요약통계량과 제시 포맷을 참고해서 필요로 하는 방법을 찾아서 사용하기 바랍니다.

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

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

728x90

저작자표시 비영리 변경금지

'R 분석과 프로그래밍 > R 통계분석' 카테고리의 다른 글

R (1) 독립된 두 표본의 모평균 차이에 대한 추정과 검정 : t.test() (11)	2015.10.15
R 두 연속형 변수 상관관계 : 공분산(covariance) cov(), 상관계수(correlation coefficient) cor() (4)	2015.10.12
R 연속형 변수 (1) 요약통계 한번에 보기 : summary(), stat.desc(), describe() (2)	2015.10.12
R 일변량 연속형 자료 요약 - (3) 분포형태와 대칭정도(distrubution, symmetry) : 왜도(skewness), 첨도(kurtosis), 분위수-분위수(Quantile-Qunantile) (0)	2015.10.12
R 일변량 연속형 자료 요약 - (2) 퍼짐 정도(dispersion) : 분산(variance), 표준편차(sd), 변이계수, 범위(range), IQR, 백분위수(percentile) (0)	2015.10.11

Posted by Rfriend

R 연속형 변수 (1) 요약통계 한번에 보기 : summary(), stat.desc(), describe()

R 분석과 프로그래밍/R 통계분석 2015. 10. 12. 16:34

지난 포스팅에서 연속형 변수(continuous variable)의 중심화 경향(central tendancy), 퍼짐 정도(dispersion), 분포형태 및 대칭정도 (distribution and symmetry) 등의 R 함수에 대해 알아보았습니다.

(1) 연속형 변수 요약통계 한번에 보기

: summary(), stat.desc(), describe()

(2) 연속형 변수 그룹별(요인별) 요약통계 비교하기

: tapply(var, factor, summary), aggregate(), summaryBy(), describe.by()

로 나누어서 소개해보겠습니다.

먼저, 한국 아침 드라마 초간단 요약 사진을 보면서 자료 정리 및 요약의 효용에 대해 직관적으로 이해를 한 후에, R을 이용한 (1) 연속형 변수 요약통계 한번에 보기를 시작하겠습니다.

[ 한국 아침 드라마 초간단 요약, ㅋㅋ~ ]

(* 사진출처 : 미상)

[ 연속형 변수 요약통계 한번에 보기 package & function 요약 ]

category	package	function	statistics
(1) descriptive statistics summary	base	summary()	min, 1Q, median, mean, 3Q, max
	pastecs	stat.desc()	nbr.val, nbr.null, nbr.na, min, max, range, sum, median, mean, var, std.dev, coef.var, skewness, kurtosis, normtest.W, normtest.p, SE.mean, CI.mean.p
	psych	describe()	n, mean, std.dev, median, trimmed, mad, min, max, range, skew, kurtosis, se
(2) descriptive statistics comparison by group (factor)	base	tapply(var, factor, summary)	summary: min, 1Q, median, mean, 3Q, max, (user defined) functions...
	base	by()	(user defined) functions...
	stats	aggregate()	(user defined) functions...
	doBy	summaryBy()	(user defined) functions...
	psych	describeBy()	n, mean, sd, median, trimmed, mad, min, max, range, skew, kurtosis, se

MASS package 내 Cars93 데이터프레임의 가격(Price), 고속도로연비(MPG.highway) 연속형 변수에 대해서 요약통계를 계산해보겠습니다.

> library(MASS)
> str(Cars93)
'data.frame':	93 obs. of  27 variables:
 $ Manufacturer      : Factor w/ 32 levels "Acura","Audi",..: 1 1 2 2 3 4 4 4 4 5 ...
 $ Model             : Factor w/ 93 levels "100","190E","240",..: 49 56 9 1 6 24 54 74 73 35 ...
 $ Type              : Factor w/ 6 levels "Compact","Large",..: 4 3 1 3 3 3 2 2 3 2 ...
 $ Min.Price         : num  12.9 29.2 25.9 30.8 23.7 14.2 19.9 22.6 26.3 33 ...
 $ Price             : num  15.9 33.9 29.1 37.7 30 15.7 20.8 23.7 26.3 34.7 ...
 $ Max.Price         : num  18.8 38.7 32.3 44.6 36.2 17.3 21.7 24.9 26.3 36.3 ...
 $ MPG.city          : int  25 18 20 19 22 22 19 16 19 16 ...
 $ MPG.highway       : int  31 25 26 26 30 31 28 25 27 25 ...
 $ AirBags           : Factor w/ 3 levels "Driver & Passenger",..: 3 1 2 1 2 2 2 2 2 2 ...
 $ DriveTrain        : Factor w/ 3 levels "4WD","Front",..: 2 2 2 2 3 2 2 3 2 2 ...
 $ Cylinders         : Factor w/ 6 levels "3","4","5","6",..: 2 4 4 4 2 2 4 4 4 5 ...
 $ EngineSize        : num  1.8 3.2 2.8 2.8 3.5 2.2 3.8 5.7 3.8 4.9 ...
 $ Horsepower        : int  140 200 172 172 208 110 170 180 170 200 ...
 $ RPM               : int  6300 5500 5500 5500 5700 5200 4800 4000 4800 4100 ...
 $ Rev.per.mile      : int  2890 2335 2280 2535 2545 2565 1570 1320 1690 1510 ...
 $ Man.trans.avail   : Factor w/ 2 levels "No","Yes": 2 2 2 2 2 1 1 1 1 1 ...
 $ Fuel.tank.capacity: num  13.2 18 16.9 21.1 21.1 16.4 18 23 18.8 18 ...
 $ Passengers        : int  5 5 5 6 4 6 6 6 5 6 ...
 $ Length            : int  177 195 180 193 186 189 200 216 198 206 ...
 $ Wheelbase         : int  102 115 102 106 109 105 111 116 108 114 ...
 $ Width             : int  68 71 67 70 69 69 74 78 73 73 ...
 $ Turn.circle       : int  37 38 37 37 39 41 42 45 41 43 ...
 $ Rear.seat.room    : num  26.5 30 28 31 27 28 30.5 30.5 26.5 35 ...
 $ Luggage.room      : int  11 15 14 17 13 16 17 21 14 18 ...
 $ Weight            : int  2705 3560 3375 3405 3640 2880 3470 4105 3495 3620 ...
 $ Origin            : Factor w/ 2 levels "USA","non-USA": 2 2 2 2 2 1 1 1 1 1 ...
 $ Make              : Factor w/ 93 levels "Acura Integra",..: 1 2 4 3 5 6 7 9 8 10 ...

각 package별 R 함수 사용예는 아래와 같습니다.

(1) base package : summary()

> ##--- (1) base package : summary()
> 
> # summary of one variable
> summary(Cars93$Price)
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
      7      12      18      20      23      62 
> 
> # summary of muliple variables
> summary(Cars93[c("Price", "MPG.highway")])
     Price     MPG.highway
 Min.   : 7   Min.   :20  
 1st Qu.:12   1st Qu.:26  
 Median :18   Median :28  
 Mean   :20   Mean   :29  
 3rd Qu.:23   3rd Qu.:31  
 Max.   :62   Max.   :50

base package의 summary() 함수는 중심화 경향과 퍼짐정도에 대해서 quick 하게 볼 수 있는 통계량들을 제공합니다. base package에 들어있으므로 사용하기에 편리한 반면, 약간 부족한 듯한 면이 있습니다.

summary() : min, 1Q, median, mean, 3Q, max

(2) pastecs package : stat.desc()

> ##--- (2) pastecs package : stat.desc()

> install.packages("pastecs")
Installing package into ‘C:/Users/user/Documents/R/win-library/3.2’
(as ‘lib’ is unspecified)
trying URL 'http://cran.rstudio.com/bin/windows/contrib/3.2/pastecs_1.3-18.zip'
Content type 'application/zip' length 1636171 bytes (1.6 MB)
downloaded 1.6 MB

package ‘pastecs’ successfully unpacked and MD5 sums checked

The downloaded binary packages are in
	C:\Users\user\AppData\Local\Temp\RtmpILyJlg\downloaded_packages
> library(pastecs)
필요한 패키지를 로딩중입니다: boot

다음의 패키지를 부착합니다: ‘boot’

The following object is masked from ‘package:survival’:

    aml

> # summary statistics of "Price", "MPG.highway" by stat.desc() function

> round(
+   stat.desc(Cars93[c("Price", "MPG.highway")], # dataframe[c(variable name)]
+                 basic = TRUE, # nbr.val, nbr.null, nbr.na, min, max, range, sum
+                 desc = TRUE, # median, mean, var, std.dev, coef.var
+                 norm = TRUE, # skewness, kurtosis, normtest.W, normtest.p
+                 p = 0.90 # SE.mean, CI.mean.p
+             )
+   , 2)
              Price MPG.highway
nbr.val       93.00       93.00
nbr.null       0.00        0.00
nbr.na         0.00        0.00
min            7.40       20.00
max           61.90       50.00
range         54.50       30.00
sum         1814.40     2705.00
median        17.70       28.00
mean          19.51       29.09
SE.mean        1.00        0.55
CI.mean.0.9    1.66        0.92
var           93.30       28.43
std.dev        9.66        5.33
coef.var       0.50        0.18
skewness       1.48        1.19
skew.2SE       2.97        2.38
kurtosis       3.05        2.30
kurt.2SE       3.08        2.32
normtest.W     0.88        0.92
normtest.p     0.00        0.00

stat.desc() 함수의 옵션 별로 제공하는 통계량은 아래와 같습니다. IQR, quantile 빼고 다 있구요, 화면에 제시해주는 포맷도 아주 깔끔해서 보기에 좋아서 제가 제일 선호하는 함수이기도 합니다.

basic = TRUE : 관측치 개수, null 개수, NA 개수, 최소값, 최대값, 범위, 합
desc = TRUE : 중앙값, 평균, 분산, 표준편차, 변이계수
norm = TRUE : 왜도, 첨도, 정규성 검정통계량, 정규성 검정 P-value
p = 0.90 : 신뢰계수 90% (유의수준 10%) 값 => 90% 신뢰구간은 평균 ± CI.mean.0.9 값
(위의 예 Price의 90% 신뢰구간은 19.51 ± 1.66)

(3) psych package : describe()

 > install.packages("psych")

> library(psych)

> describe(Cars93[c("Price", "MPG.highway")], 
+          na.rm = TRUE, # missing value is not included 
+          interp = TRUE, # method of median calculation
+          skew = TRUE, # skewness, kurtosis
+          ranges = TRUE, # ragne
+          trim = 0.1 # trimmed mean
+          )
           vars n mean  sd median trimmed  mad  min  max range skew kurtosis se
Price       1  93   20 9.7     18      18   8.3  7.4  62   54   1.5   3.0  1.00
MPG.highway 2  93   29 5.3     28      29   4.4 20.0  50   30   1.2   2.3  0.55

psych package의 describe() 함수는 summary() 보다는 많고 stat.desc() 보다는 적은 기술통계량을 보여줍니다.

describe() : 관측값 개수(n), 평균(mean), 표준편차(sd), 중앙값(median), 절삭평균(10% 절삭평균), 중위값절대편차(from 중위값) (MAD, median absolute deviation), 최소값(min), 최대값(max), 범위(range), 왜도(skew), 첨도(kurtosis), 표준오차(SE, standard error)

[참고]

중위값 절대 편차 (MAD, median absolute deviation) = median(|X - median(X)|) * 1.4826
(이때 1.4826은 scaling factor (또는 normalizing constant) 이며, 정규적인 자료에서 scaling factor를 곱해주면 표준편차와 비슷해짐)

(* 댓글남겨주신 데이터과학자님, 감사드립니다)

summary()와 stat.desc() 함수가 column에 variable을 row에 statistics를 보여주는 반면, describe() 함수는 반대로 row에 variable을, column에 statistics를 보여줍니다.

원하는 제시 포맷을 선택하고, 원하는 기술통계량을 옵션에서 선택해서 연속형 자료를 일목요연하게 정리, 요약하면 되겠습니다.

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

질문은 댓글로 남겨주세요.

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

728x90

저작자표시 비영리 변경금지

'R 분석과 프로그래밍 > R 통계분석' 카테고리의 다른 글

R 두 연속형 변수 상관관계 : 공분산(covariance) cov(), 상관계수(correlation coefficient) cor() (4)	2015.10.12
R 연속형 변수 (2) 그룹별(요인별) 요약통계 비교하기 : tapply(var, factor, summary), by(), aggregate(), summaryBy(), describeBy() (4)	2015.10.12
R 일변량 연속형 자료 요약 - (3) 분포형태와 대칭정도(distrubution, symmetry) : 왜도(skewness), 첨도(kurtosis), 분위수-분위수(Quantile-Qunantile) (0)	2015.10.12
R 일변량 연속형 자료 요약 - (2) 퍼짐 정도(dispersion) : 분산(variance), 표준편차(sd), 변이계수, 범위(range), IQR, 백분위수(percentile) (0)	2015.10.11
R 일변량 연속형 자료 요약 - (1) 중심화 경향 : 산술평균(mean), 중앙값(median), 최빈값, 기하평균, CAGR, 조화평균, 가중평균 (2)	2015.10.11

Posted by Rfriend

R 일변량 연속형 자료 요약 - (3) 분포형태와 대칭정도(distrubution, symmetry) : 왜도(skewness), 첨도(kurtosis), 분위수-분위수(Quantile-Qunantile)

R 분석과 프로그래밍/R 통계분석 2015. 10. 12. 10:38

일변량 연속형 자료에 대해 기술통계량(descriptive statistics)을 이용한 자료의 요약과 정리는 크게

- (1) 중심화 경향 (central tendency)

: 산술평균, 중앙값, 최빈값, 기하평균, CAGR, 조화평균, 가중평균

- (2) 퍼짐 정도 (dispersion)

: 분산, 표준편차, 변이계수, 범위, IQR, 백분위수

- (3) 분포형태와 대칭정도 (distribution)

: 왜도, 첨도, 분위수-분위수

의 3가지로 구분할 수 있습니다.

이번 포스팅에서는 일변량 연속형 자료의 (3) 분포형태와 대칭 정도 (distribution and symmetry)에 대해 통계 이론과 활용 상의 주의점을 알아보고, R 함수를 가지고 예를 들어보겠습니다.

다시 한번 말씀드리지만, 위의 3개의 카테고리는 설명을 쉽게 하기 위해 구조화한 것이며, 실제 분석업무할 때는 중심 경향, 퍼짐 정도, 분포형태와 대칭 정도의 통계량을 같이 봐야 하며, 그래프도 병행해서 보면서 해석을 하는 것이 중요합니다. 통계량을 단편적으로 보거나, 통계량만 보는 것은 잘못된 해석, 판단으로 이끌 수 있는 위험이 있습니다.

[ 산술통계량(descriptive statistics)과 R function ]

산술통계	통계량 (statistics)	R function
중심화 경향 (central tendency)	산술평균 (arithmetic mean)	mean()
	중앙값 (median)	median()
	최빈값 (mode)	which.max(table())
	기하평균 (geometric mean)	prod(x)^(1/n)1/mean(1/x) where, n = length(x)
	연평균성장률 (CAGR : Componded Average Growth Rate)	(FV/IV)^(1/n)-1 where, IV : initial value of an investment FV : final value of an investment n : investment periods
	조화평균 (harmonic mean)	1/mean(1/x)
	가중평균 (weighted average)	weighted.mean()


퍼짐 정도 (dispersion)	분산 (variance)	var()
	표준편차 (standard deviation)	sd()
	변이계수 (coefficient of variation)	100*sd(x)/mean(x)
	범위 (range)	diff(range())
	IQR (Inter Quartile Range)	IQR()
	최소값 (min)	min()
	최대값 (max)	max()
	백분위수(percentile)	quantile(x, probs=c(,,,,))
분포형태와 대칭정도 (distribution, symmetry)	왜도 (skewness)	skewness(), fBasics package
	첨도 (kurtosis)	kurtosis(), fBasics package
	분위수-분위수(Quantile-Quantile)	qqnorm(), qqline(), qqplot()

※ 중심화 경향, 퍼짐 정도, 분포형태와 대칭정도의 통계량을 함께 봐야함

※ 통계량과 함께 그래프를 함께 봐야함

R 실습에는 MASS 패키지 내 Cars93 데이터의 차종(Type), 가격(Price) 변수를 활용하겠습니다.

> library(MASS)

> str(Cars93)
'data.frame':	93 obs. of  27 variables:
 $ Manufacturer      : Factor w/ 32 levels "Acura","Audi",..: 1 1 2 2 3 4 4 4 4 5 ...
 $ Model             : Factor w/ 93 levels "100","190E","240",..: 49 56 9 1 6 24 54 74 73 35 ...
 $ Type              : Factor w/ 6 levels "Compact","Large",..: 4 3 1 3 3 3 2 2 3 2 ...
 $ Min.Price         : num  12.9 29.2 25.9 30.8 23.7 14.2 19.9 22.6 26.3 33 ...
 $ Price             : num  15.9 33.9 29.1 37.7 30 15.7 20.8 23.7 26.3 34.7 ...
 $ Max.Price         : num  18.8 38.7 32.3 44.6 36.2 17.3 21.7 24.9 26.3 36.3 ...
 $ MPG.city          : int  25 18 20 19 22 22 19 16 19 16 ...
 $ MPG.highway       : int  31 25 26 26 30 31 28 25 27 25 ...
 $ AirBags           : Factor w/ 3 levels "Driver & Passenger",..: 3 1 2 1 2 2 2 2 2 2 ...
 $ DriveTrain        : Factor w/ 3 levels "4WD","Front",..: 2 2 2 2 3 2 2 3 2 2 ...
 $ Cylinders         : Factor w/ 6 levels "3","4","5","6",..: 2 4 4 4 2 2 4 4 4 5 ...
 $ EngineSize        : num  1.8 3.2 2.8 2.8 3.5 2.2 3.8 5.7 3.8 4.9 ...
 $ Horsepower        : int  140 200 172 172 208 110 170 180 170 200 ...
 $ RPM               : int  6300 5500 5500 5500 5700 5200 4800 4000 4800 4100 ...
 $ Rev.per.mile      : int  2890 2335 2280 2535 2545 2565 1570 1320 1690 1510 ...
 $ Man.trans.avail   : Factor w/ 2 levels "No","Yes": 2 2 2 2 2 1 1 1 1 1 ...
 $ Fuel.tank.capacity: num  13.2 18 16.9 21.1 21.1 16.4 18 23 18.8 18 ...
 $ Passengers        : int  5 5 5 6 4 6 6 6 5 6 ...
 $ Length            : int  177 195 180 193 186 189 200 216 198 206 ...
 $ Wheelbase         : int  102 115 102 106 109 105 111 116 108 114 ...
 $ Width             : int  68 71 67 70 69 69 74 78 73 73 ...
 $ Turn.circle       : int  37 38 37 37 39 41 42 45 41 43 ...
 $ Rear.seat.room    : num  26.5 30 28 31 27 28 30.5 30.5 26.5 35 ...
 $ Luggage.room      : int  11 15 14 17 13 16 17 21 14 18 ...
 $ Weight            : int  2705 3560 3375 3405 3640 2880 3470 4105 3495 3620 ...
 $ Origin            : Factor w/ 2 levels "USA","non-USA": 2 2 2 2 2 1 1 1 1 1 ...
 $ Make              : Factor w/ 93 levels "Acura Integra",..: 1 2 4 3 5 6 7 9 8 10 .

Histogram으로 살펴본 차종별 가격(Price by Type)의 분포는 아래와 같으며, 그래프를 통해 차종별 가격의 왜도(좌우 대칭 정도)와 첨도(정규분포 대비 봉우리 높이 정도)를 가늠해 볼 수 있습니다.

> # Histogram, Price by Car Type

> library(ggplot2)
> ggplot(Cars93, aes(x=Price)) + 
+   geom_histogram(binwidth=3, fill = "blue", colour = "black") + 
+   ggtitle("Histogram, Price by Type") + 
+   facet_grid(Type ~ .)

fBasics package를 먼저 설치, 호출한 후에 왜도 skewness() 함수, 첨도 kurtosis() 함수를 사용해보겠습니다.

> install.packages("fBasics")
Installing package into ‘C:/Users/user/Documents/R/win-library/3.2’
(as ‘lib’ is unspecified)
trying URL 'http://cran.rstudio.com/bin/windows/contrib/3.2/fBasics_3011.87.zip'
Content type 'application/zip' length 1553624 bytes (1.5 MB)
downloaded 1.5 MB

package ‘fBasics’ successfully unpacked and MD5 sums checked

The downloaded binary packages are in
	C:\Users\user\AppData\Local\Temp\Rtmpu6UbP8\downloaded_packages
> library(fBasics)
필요한 패키지를 로딩중입니다: timeDate
필요한 패키지를 로딩중입니다: timeSeries


Rmetrics Package fBasics
Analysing Markets and calculating Basic Statistics
Copyright (C) 2005-2014 Rmetrics Association Zurich
Educational Software for Financial Engineering and Computational Science
Rmetrics is free software and comes with ABSOLUTELY NO WARRANTY.
https://www.rmetrics.org --- Mail to: info@rmetrics.org

R의 fBasics 패키지를 활용해 왜도(skewness), 첨도(kurtosis)를 분석해보겠습니다.

(1) 왜도 (skewness) : skewness()

왜도(skewness)는 자료의 대칭성을 알아보는 측도입니다. 좌우 대칭형태를 띠는 정규분포는 왜도(β3) 점수가 '0' 이며, 오른쪽으로 꼬리가 긴 분포 (right-skewed distribution)은 왜도 점수가 '0'보다 크며(β3 > 0), 왼쪽으로 꼬리가 긴 분포(left-skewed distribution)은 왜도 점수가 '0'보다 작습니다((β3 < 0). 이는 위의 수식을 보면 기본적으로 관측값에서 평균을 뺀 값을 세제곱한 후에 더했기 때문입니다.

위의 Histogram을 보면 Cars93의 Type별 Price를 보면 전반적으로 오른쪽으로 꼬리가 긴 형태를 띠고 있는데요, 아래의 왜도를 구한 값을 보면 모두 양수 임을 알 수 있고, 특히 오른쪽으로 꼬리가 긴 Midsize와 Sporty type은 왜도 점수가 높음을 알 수 있습니다.

> # skewness : skewness()
> skewness(Cars93$Price)
[1] 1.483982
attr(,"method")
[1] "moment"
>

> with(Cars93, tapply(Price, Type, skewness))
  Compact     Large   Midsize     Small    Sporty       Van 
0.7797609 0.8152686 1.0600340 1.0068697 1.0286200 0.1655290

(2) 첨도 (kurtosis) : kurtosis()

첨도(kurtosis)는 정규분포 대비 봉오리의 높이를 알아보는 측도입니다. 첨도가 '0'보다 크면 (β4 > 0) 정규분

포보다 뾰족하다는 의미이며, 첨도가 '0'보다 작으면 (β4 < 0) 정규분포보다 납작하다는 뜻으로 해석하면 되겠습니다.

더불어, 첨도(kurtosis, β4)는 데이터가 이봉분포(two mountaintop distribution)에 대해서 얼마나 단봉분포(one mountaintop distribution)에 가깝게 있는가를 판단하는데도 사용합니다.

아래에는 차종(Type)별 가격(Price)의 첨도를 계산해 놓았는데요, Midsize와 Small type이 양수로 나와서 정규분포보다 뾰족한 형태를 취하고 있다고 볼 수 있겠습니다. (위의 histogram 을 봤을 때, Van이 음수로 나온게 의외네요)

> # kurtosis : kurtosis()
> 
> kurtosis(Cars93$Price)
[1] 3.051418
attr(,"method")
[1] "excess"
> 
> with(Cars93, tapply(Price, Type, kurtosis))
   Compact      Large    Midsize      Small     Sporty        Van 
-0.8006747 -1.0263810  0.6455252  1.0822979 -0.1019349 -0.5974974

마지막으로 정규분포 형태를 띠는지의 검정 및 그래프 확인하는 방법은 아래의 포스팅을 참조하시기 바랍니다.

☞ R 단일 모집단 분포의 정규성 검정 : shapiro.test(), qqnorm(), qqline()

☞ R 데이터 변환 (2) 정규분포화 log(), sqrt()

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

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

728x90

저작자표시 비영리 변경금지

'R 분석과 프로그래밍 > R 통계분석' 카테고리의 다른 글

R 연속형 변수 (2) 그룹별(요인별) 요약통계 비교하기 : tapply(var, factor, summary), by(), aggregate(), summaryBy(), describeBy() (4)	2015.10.12
R 연속형 변수 (1) 요약통계 한번에 보기 : summary(), stat.desc(), describe() (2)	2015.10.12
R 일변량 연속형 자료 요약 - (2) 퍼짐 정도(dispersion) : 분산(variance), 표준편차(sd), 변이계수, 범위(range), IQR, 백분위수(percentile) (0)	2015.10.11
R 일변량 연속형 자료 요약 - (1) 중심화 경향 : 산술평균(mean), 중앙값(median), 최빈값, 기하평균, CAGR, 조화평균, 가중평균 (2)	2015.10.11
R 이변량 범주형 자료의 정리와 요약 - 분할표(contingency table) : table(), xtabs(), prop.table(), margin.table(), addmargins(), CrossTable() (0)	2015.10.10

Posted by Rfriend

R 일변량 연속형 자료 요약 - (2) 퍼짐 정도(dispersion) : 분산(variance), 표준편차(sd), 변이계수, 범위(range), IQR, 백분위수(percentile)

R 분석과 프로그래밍/R 통계분석 2015. 10. 11. 23:30

일변량 연속형 자료에 대해 기술통계량(descriptive statistics)을 이용한 자료의 요약과 정리는 크게

- (1) 중심화 경향 (central tendency)

: 산술평균, 중앙값, 최빈값, 기하평균, CAGR, 조화평균, 가중평균

- (2) 퍼짐 정도 (dispersion)

: 분산, 표준편차, 변이계수, 범위, IQR, 백분위수

- (3) 분포형태와 대칭정도 (distribution)

: 왜도, 첨도, 분위수-분위수

의 3가지로 구분할 수 있습니다.

지난 포스팅에서는 중심화 경향에 대해서 알아보았는데요, 이것만 가지고는 자료의 특성을 파악했다고 보기 어려우며, 이와 더불어 자료가 중심으로 부터 얼마나 퍼져있는지, 분포는 어떤 형태인지를 같이 알아야만 합니다.

아래 3-1반과 3-2반의 수학 점수를 보면 두 학급 모두 평균은 62점으로 같습니다만, 표준편차는 27점 vs. 5.7점으로 매우 다름을 알 수 있습니다. 3-1반은 최우등생과 최열등생이 모여있는 반이고, 3-2반은 비슷한 실력의 중급 학생들이 모여있는 반이라고 하겠습니다. 왜 평균만 보면 안되는지 아셨을 겁니다.

학급 (class)	수학 점수 (math score)	평균 (mean)	표준편차(sd)
3학년 1반	25, 55, 60, 70, 100	62	27.06
3학년 2반	55, 60, 60, 65, 70	62	5.70

이번 포스팅에서는 일변량 연속형 자료의 (2) 퍼짐 정도 (dispersion)에 대해 통계 이론과 활용 상의 주의점을 알아보고, R 함수를 가지고 예를 들어보겠습니다.

[ 산술통계량(descriptive statistics)과 R function ]

산술통계	통계량 (statistics)	R function
중심화 경향 (central tendency)	산술평균 (arithmetic mean)	mean()
	중앙값 (median)	median()
	최빈값 (mode)	which.max(table())
	기하평균 (geometric mean)	prod(x)^(1/n)1/mean(1/x) where, n = length(x)
	연평균성장률 (CAGR : Componded Average Growth Rate)	(FV/IV)^(1/n)-1 where, IV : initial value of an investment FV : final value of an investment n : investment periods
	조화평균 (harmonic mean)	1/mean(1/x)
	가중평균 (weighted average)	weighted.mean()


퍼짐 정도 (dispersion)	분산 (variance)	var()
	표준편차 (standard deviation)	sd()
	변이계수 (coefficient of variation)	*100sd(x)/mean(x)**
	범위 (range)	diff(range())
	IQR (Inter Quartile Range)	IQR()
	최소값 (min)	min()
	최대값 (max)	max()
	백분위수(percentile)	quantile(x, probs=c(,,,,))
분포형태와 대칭정도 (distribution)	왜도 (skewness)	skewness(), fBasics package
	첨도 (kurtosis)	kurtosis(), fBasics package
	분위수-분위수(Quantile-Quantile)	qqnorm(), qqline(), qqplot()

※ 중심화 경향, 퍼짐 정도, 분포형태와 대칭정도의 통계량을 함께 봐야함

※ 통계량과 함께 그래프를 함께 봐야함

R 실습에는 MASS 패키지 내 Cars93 데이터의 차종(Type), 가격(Price) 변수를 활용하겠습니다.

> library(MASS)

> str(Cars93)
'data.frame':	93 obs. of  27 variables:
 $ Manufacturer      : Factor w/ 32 levels "Acura","Audi",..: 1 1 2 2 3 4 4 4 4 5 ...
 $ Model             : Factor w/ 93 levels "100","190E","240",..: 49 56 9 1 6 24 54 74 73 35 ...
 $ Type              : Factor w/ 6 levels "Compact","Large",..: 4 3 1 3 3 3 2 2 3 2 ...
 $ Min.Price         : num  12.9 29.2 25.9 30.8 23.7 14.2 19.9 22.6 26.3 33 ...
 $ Price             : num  15.9 33.9 29.1 37.7 30 15.7 20.8 23.7 26.3 34.7 ...
 $ Max.Price         : num  18.8 38.7 32.3 44.6 36.2 17.3 21.7 24.9 26.3 36.3 ...
 $ MPG.city          : int  25 18 20 19 22 22 19 16 19 16 ...
 $ MPG.highway       : int  31 25 26 26 30 31 28 25 27 25 ...
 $ AirBags           : Factor w/ 3 levels "Driver & Passenger",..: 3 1 2 1 2 2 2 2 2 2 ...
 $ DriveTrain        : Factor w/ 3 levels "4WD","Front",..: 2 2 2 2 3 2 2 3 2 2 ...
 $ Cylinders         : Factor w/ 6 levels "3","4","5","6",..: 2 4 4 4 2 2 4 4 4 5 ...
 $ EngineSize        : num  1.8 3.2 2.8 2.8 3.5 2.2 3.8 5.7 3.8 4.9 ...
 $ Horsepower        : int  140 200 172 172 208 110 170 180 170 200 ...
 $ RPM               : int  6300 5500 5500 5500 5700 5200 4800 4000 4800 4100 ...
 $ Rev.per.mile      : int  2890 2335 2280 2535 2545 2565 1570 1320 1690 1510 ...
 $ Man.trans.avail   : Factor w/ 2 levels "No","Yes": 2 2 2 2 2 1 1 1 1 1 ...
 $ Fuel.tank.capacity: num  13.2 18 16.9 21.1 21.1 16.4 18 23 18.8 18 ...
 $ Passengers        : int  5 5 5 6 4 6 6 6 5 6 ...
 $ Length            : int  177 195 180 193 186 189 200 216 198 206 ...
 $ Wheelbase         : int  102 115 102 106 109 105 111 116 108 114 ...
 $ Width             : int  68 71 67 70 69 69 74 78 73 73 ...
 $ Turn.circle       : int  37 38 37 37 39 41 42 45 41 43 ...
 $ Rear.seat.room    : num  26.5 30 28 31 27 28 30.5 30.5 26.5 35 ...
 $ Luggage.room      : int  11 15 14 17 13 16 17 21 14 18 ...
 $ Weight            : int  2705 3560 3375 3405 3640 2880 3470 4105 3495 3620 ...
 $ Origin            : Factor w/ 2 levels "USA","non-USA": 2 2 2 2 2 1 1 1 1 1 ...
 $ Make              : Factor w/ 93 levels "Acura Integra",..: 1 2 4 3 5 6 7 9 8 10 .

아래의 차종별로 가격 (Price by Type)을 Histogram으로 살펴보면 Midsize 가 좌우로 가장 많이 퍼져있으며, Compact, Large, Sporty 가 그 다음으로 많이 퍼져있고, 마지막으로 Small, Van 이 가장 작게 퍼져있음을 알 수 있습니다. 이걸 아래의 통계량들을 가지고 퍼짐 정도를 측정해 보겠습니다.

> # Histogram, Price by Type
> library(MASS) # Cars93 dataset
> library(ggplot2)
> ggplot(Cars93, aes(x=Price)) + 
+   geom_histogram(binwidth=5, fill = "blue", colour = "black") + 
+   ggtitle("Histogram, Price by Type") + 
+   facet_grid(Type ~ .)

이제부터 R 함수를 이용해서 퍼짐 정도 (dispersion) 를 파악할 수 있는 통계량을 하나씩 살펴보겠습니다.

(1) 분산 (variance) : var()

분산(variance)은 표준편차(standard deviation)와 함께 가장 일반적으로 사용되는 퍼짐 정도를 나타내는 통계량입니다. 각 관찰값에서 평균을 빼면 평균으로 부터의 거리, 편차(deviation)가 나오는데요, 이걸 모두 합하면 '0'이 됩니다. 따라서 '0'이 되지 않고 퍼진 정도를 알기 위해서 제곱(square)을 하여 합(sum)한 것이고, 관찰값 개수 N으로 나누어서 편차제곱의 평균값으로 퍼진 정도를 측정한 것이 분산(variance)입니다.

표본에서 분산을 계산할 때는 편차 제곱합을 관찰값 개수 n에서 1을 뺀 n-1을 사용하여 나누어줍니다.

> # variance : var()
> 
> var(Cars93$Price)
[1] 93.30458
> 
> with(Cars93, tapply(Price, Type, var))
   Compact      Large    Midsize      Small     Sporty        Van 
 44.714500  40.164000 150.426320   3.815333  63.596099   3.527500

차종별 가격(Price by Type)의 분산을 구하기 위해 tapply(var, factor, function) 함수를 사용하였습니다.

(2) 표준편차 (standard deviation) : sd()

표준편차(standard deviation)는 분산(variance)에다가 제곱근(squared root)을 취한 값입니다. 분산(variance)의 경우 편차를 제곱하다 보니 원자료의 scale과는 달라져버리게 되어 해석하는데 좀 곤란한 상황이 벌어집니다. 이 문제를 해결할 수 있는 것이 바로 표준편차입니다. 편차 제곱한 분산에다가 제곱근을 취했기 때문에 원자료와 scale이 동일해지기 때문입니다. 표준편차도 분산과 동일하게 숫자가 커질 수록 중심으로부터 멀리 퍼져있다고 해석하면 되며, 원자료와 scale이 동일하기 때문에 평균에서 (정규분포의 경우) 좌우로 표준편차만큼 퍼져있다고 생각하면 이해하기가 쉽겠습니다.

> # standard deviation : sd()
> 
> sd(Cars93$Price)
[1] 9.65943
> 
> with(Cars93, tapply(Price, Type, sd))
  Compact     Large   Midsize     Small    Sporty       Van 
 6.686890  6.337507 12.264841  1.953288  7.974716  1.878164

위의 차종별 가격의 표준편차를 보면 위의 histogram과 동일한 결과가 나왔음을 알 수 있습니다. Midsize가 표준편차가 12.26으로 가장 크고, Van이 1.87로 표준편차가 가장 작게 나왔습니다.

(3) 변이계수 (coefficeint of variation) : 100*sd()/mean()

위에서 표준편차(standard deviation)가 scale이 원자료와 같기 때문에 분산(variance)보다는 사용하기에 유용하다고 말했습니다. 하지만 표준편차도 약점이 있는데요, 절대 크기가 현저하게 달라서 평균이 서로 매우 다른 두 집단 간 비교, 측정 단위가 다른 두 변수 간 비교에는 부적합합니다. 이럴 때 퍼짐 정도를 비교 가능하도록 표준화해준 통계량이 변이계수(coeffieicent of variation)이 되겠습니다. 변이계수는 표준편차를 평균으로 나눈 다음에 100을 곱해서 계산합니다.

차종별 가격의 변이계수를 구하면 아래와 같은데요, 변이계수가 표준편차와 뭐가 다른가 잘 감이 안잡힐 수도 있겠습니다.

> # coefficient of variation : sd()/mean()
> 
> with(Cars93, 100*sd(Price)/mean(Price))
[1] 49.51096
> 
> attach(Cars93)
> with(Cars93[Type == c("Compact"),], 100*sd(Price)/mean(Price))
[1] 36.71594
> with(Cars93[Type == c("Large"),], 100*sd(Price)/mean(Price))
[1] 26.08028
> with(Cars93[Type == c("Midsize"),], 100*sd(Price)/mean(Price))
[1] 45.06121
> with(Cars93[Type == c("Small"),], 100*sd(Price)/mean(Price))
[1] 19.21267
> with(Cars93[Type == c("Sporty"),], 100*sd(Price)/mean(Price))
[1] 41.12193
> with(Cars93[Type == c("Van"),], 100*sd(Price)/mean(Price))
[1] 9.833319
> detach(Cars93)

변이계수의 이해를 돕기 위해서 하나의 예를 추가로 들어보겠습니다.

A회사와 B회사가 있는데요, 한달 주식가격의 평균과 표준편차가 아래와 같은 때, ~~표준편차로만 보면 B회사(sd 2,000원)가 A회사(sd 1,000원)의 2배로서 Risk가 더 높다고 생각할 수 있습니다만~~, 여기에는 함정이 있으며, 이렇게 계산하면 틀립니다. B회사의 주당 평균 주가(mean 50,000원)는 A회사의 주당 평균주가(mean 10,000원)의 5배에 해당할만큼 큰 차이를 보이고 있습니다.

이럴 경우 급이 다르기 때문에 평균으로 표준편차를 나누어준 비율인 변이계수를 사용해서 동급으로 만들어주고 퍼짐 정도를 비교해야만 합니다. A회사의 변이계수는 10%, B회사의 변이계수는 4%로서 A회사가 B회사보다 Risk가 2.5배 더 높다고 평가할 수 있으며, 앞서의 표준편차와는 정반대의 결과가 나왔음에 유의하시기 바랍니다.

> # example : stock price's mean, sd of company A and company B
> 
> company_A_mean <- c(10000)
> company_A_sd <- c(1000)
> 
> company_B_mean <- c(50000)
> company_B_sd <- c(2000)
> 
> 
> coe_var_A <- 100*company_A_sd/company_A_mean
> coe_var_A
[1] 10
> 
> coe_var_B <- 100*company_B_sd/company_B_mean
> coe_var_B
[1] 4

(4) 최소값 (min) : min()

(5) 최대값 (max) : max()

(6) 범위 (range) : diff(range())

(7) 백분위수 (percentile) : quantile(x, probs=c(,,,,))

(8) IQR (Inter Quartile Range) : IQR()

범위(range)는 최대값에서 최소값을 뺀 값으로, 직관적으로 가장 이해하기 쉬운 퍼짐 정도 통계량입니다. 다만, 특이값(outlier)에 민감하므로 특이값을 제거 후에 사용하거나, 아니면 특이값에 견고한 IQR(Inter Quartile Range) 를 대신 사용할 수 있습니다.

제 p 백분위수(pth percentile)는 자료를 크기 순서대로 정렬해놓았을 때 p%가 자기값 이하(자기값 포함)로 적어도 p%의 관측값이 있고, 자기값 이상으로 적오도 (100-p)%의 관측값이 있는 수를 의미합니다. Q1, Q2(median), Q3 등은 우리가 자주 사용하는 대표적인 백분위수(percentile)로서, 사분위수(quartile)이라고도 하며 이때 Q1은 25% percentile, Q2는 50% percentile, Q3는 75% percentile이 되겠지요.

R로는 함수 한줄로 누워서 떡먹기보다 더 쉬운데요, 이것을 SQL, Hive로 구현하려면 머리가 좀 아프고 코딩을 좀 해야만 합니다. ^^;

자, 그럼 R로 차종별 가격의 Min, Max, 범위, 25% percentile(Q1), 75% percentile(Q3), IQR을 차례대로 구해보겠습니다.

> ##----------
> # min, max, range, IQR, percentile
> attach(Cars93)
> 
> # min : min()
> min(Price)
[1] 7.4
> tapply(Price, Type, min)
Compact   Large Midsize   Small  Sporty     Van 
   11.1    18.4    13.9     7.4    10.0    16.3 
> 
> # max : max()
> max(Price)
[1] 61.9
> tapply(Price, Type, max)
Compact   Large Midsize   Small  Sporty     Van 
   31.9    36.1    61.9    15.9    38.0    22.7 
> 
> # range : diff(range())
> diff(range(Price))
[1] 54.5
> 
> diff(range(Cars93[Type==c("Compact"),]$Price))
[1] 20.8
> diff(range(Cars93[Type==c("Large"),]$Price))
[1] 17.7
> diff(range(Cars93[Type==c("Midsize"),]$Price))
[1] 48
> diff(range(Cars93[Type==c("Small"),]$Price))
[1] 8.5
> diff(range(Cars93[Type==c("Sporty"),]$Price))
[1] 28
> diff(range(Cars93[Type==c("Van"),]$Price))
[1] 6.4
> 
> # Percentile : quantile(var, probs=c(,,))
> quantile(Price, c(0.25, 0.75))
 25%  75% 
12.2 23.3 
> 
> quantile(Cars93[Type==c("Compact"),]$Price, c(0.25, 0.75))
   25%    75% 
13.375 20.675 
> quantile(Cars93[Type==c("Large"),]$Price, c(0.25, 0.75))
  25%   75% 
20.00 26.95 
> quantile(Cars93[Type==c("Midsize"),]$Price, c(0.25, 0.75))
   25%    75% 
16.775 34.200 
> quantile(Cars93[Type==c("Small"),]$Price, c(0.25, 0.75))
 25%  75% 
 8.6 11.3 
> quantile(Cars93[Type==c("Sporty"),]$Price, c(0.25, 0.75))
   25%    75% 
14.175 22.425 
> quantile(Cars93[Type==c("Van"),]$Price, c(0.25, 0.75))
 25%  75% 
19.0 19.7 
> 
> 
> # IQR : IQR()
> IQR(Price)
[1] 11.1
> 
> IQR(Cars93[Type==c("Compact"),]$Price)
[1] 7.3
> IQR(Cars93[Type==c("Large"),]$Price)
[1] 6.95
> IQR(Cars93[Type==c("Midsize"),]$Price)
[1] 17.425
> IQR(Cars93[Type==c("Small"),]$Price)
[1] 2.7
> IQR(Cars93[Type==c("Sporty"),]$Price)
[1] 8.25
> IQR(Cars93[Type==c("Van"),]$Price)
[1] 0.7
> detach(Cars93)

위의 퍼짐 정도(range, Q1, median, Q3, lower/upper whisker line, outlier) & 중심 경향(mean) 관련 통계량들을 박스 그림(box-and-whis ker plot)으로 그리면 아래와 같습니다.

> # box plot with mean
> ggplot(Cars93, aes(x = Type, y = Price)) +
+ geom_boxplot(width=0.8, outlier.size=3, outlier.shape=16, outlier.colour="red") +
+ stat_summary(fun.y="mean", geom="point", shape=21, size=3, fill="blue") +
+ ggtitle("Box Plot by Car Type, adding mean")

저 위에도 적어놨지만요, 통계량은 중심화 경향, 퍼짐 정도, 분포형태 및 대칭 정도 통계량을 같이 봐야 하고, 그래프도 같이 봐서 종합적으로 해석하는 것이 정말 중요합니다.

중심화 경향과 퍼짐 정도가 다른 두 데이터셋을 표준화하는 방법은 아래의 포스팅을 참고하시기 바랍니다.

☞ R 데이터 변환 (1) 표준화 : z 표준화 변환, [0-1] 변환

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

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

728x90

저작자표시 비영리 변경금지

'R 분석과 프로그래밍 > R 통계분석' 카테고리의 다른 글

R 연속형 변수 (1) 요약통계 한번에 보기 : summary(), stat.desc(), describe() (2)	2015.10.12
R 일변량 연속형 자료 요약 - (3) 분포형태와 대칭정도(distrubution, symmetry) : 왜도(skewness), 첨도(kurtosis), 분위수-분위수(Quantile-Qunantile) (0)	2015.10.12
R 일변량 연속형 자료 요약 - (1) 중심화 경향 : 산술평균(mean), 중앙값(median), 최빈값, 기하평균, CAGR, 조화평균, 가중평균 (2)	2015.10.11
R 이변량 범주형 자료의 정리와 요약 - 분할표(contingency table) : table(), xtabs(), prop.table(), margin.table(), addmargins(), CrossTable() (0)	2015.10.10
R 일변량 범주형 자료의 정리와 요약 - 도수분포표(frequency distribution table) : table(), xtabs() (2)	2015.10.10

Posted by Rfriend

R 일변량 연속형 자료 요약 - (1) 중심화 경향 : 산술평균(mean), 중앙값(median), 최빈값, 기하평균, CAGR, 조화평균, 가중평균

R 분석과 프로그래밍/R 통계분석 2015. 10. 11. 01:16

일변량 연속형 자료에 대해 기술통계량(descriptive statistics)을 이용한 자료의 요약과 정리는 크게

- (1) 중심화 경향 (central tendency)

: 산술평균, 중앙값, 최빈값, 기하평균, CAGR, 조화평균, 가중평균

- (2) 퍼짐 정도 (dispersion)

: 분산, 표준편차, 변이계수, 범위, IQR, 백분위수

- (3) 분포형태와 대칭정도 (distribution)

: 왜도, 첨도, 분위수-분위수

의 3가지로 구분할 수 있습니다.

이번 포스팅에서는 일변량 연속형 자료의 (1) 중심화 경향에 대해 통계 이론과 활용 상의 주의점을 알아보고, R 함수를 가지고 예를 들어보겠습니다.

일반적으로 많이 사용되는 산술평균, 중앙값, 최빈값을 먼저 살펴보겠습니다.

이어서 산술평균의 함정의 주의사항과 함께 기하평균, 연평균성장률, 조화평균, 가중평균을 언제 사용해야 하는지, 어떻게 계산하는지를 차례로 알아보겠습니다.

[ 산술통계량(descriptive statistics)과 R function ]

산술통계	통계량 (statistics)	R function
중심화 경향 (central tendency)	산술평균 (arithmetic mean)	mean()
	중앙값 (median)	median()
	최빈값 (mode)	which.max(table())
	기하평균 (geometric mean)	prod(x)^(1/n)1/mean(1/x) where, n = length(x)
	연평균성장률 (CAGR : Componded Average Growth Rate)	(FV/IV)^(1/n)-1 where, IV : initial value of an investment FV : final value of an investment n : investment periods
	조화평균 (harmonic mean)	1/mean(1/x)
	가중평균 (weighted average)	weighted.mean()


퍼짐 정도 (dispersion)	분산 (variance)	var()
	표준편차 (standard deviation)	sd()
	변이계수 (coefficient of variation)	100*sd(x)/mean(x)
	범위 (range)	diff(range())
	IQR (Inter Quartile Range)	IQR()
	최소값 (min)	min()
	최대값 (max)	max()
	백분위수(percentile)	quantile(x, probs=c(,,,,))
분포형태와 대칭정도 (distribution)	왜도 (skewness)	skewness(), fBasics package
	첨도 (kurtosis)	kurtosis(), fBasics package
	분위수-분위수(Quantile-Quantile)	qqnorm(), qqline(), qqplot()

※ 중심화 경향, 퍼짐 정도, 분포형태와 대칭정도의 통계량을 함께 봐야함

※ 통계량과 함께 그래프를 함께 봐야함

R 실습에는 MASS 패키지 내 Cars93 데이터의 가격(Price) 변수를 활용하겠습니다.

> library(MASS)

> str(Cars93)
'data.frame':	93 obs. of  27 variables:
 $ Manufacturer      : Factor w/ 32 levels "Acura","Audi",..: 1 1 2 2 3 4 4 4 4 5 ...
 $ Model             : Factor w/ 93 levels "100","190E","240",..: 49 56 9 1 6 24 54 74 73 35 ...
 $ Type              : Factor w/ 6 levels "Compact","Large",..: 4 3 1 3 3 3 2 2 3 2 ...
 $ Min.Price         : num  12.9 29.2 25.9 30.8 23.7 14.2 19.9 22.6 26.3 33 ...
 $ Price             : num  15.9 33.9 29.1 37.7 30 15.7 20.8 23.7 26.3 34.7 ...
 $ Max.Price         : num  18.8 38.7 32.3 44.6 36.2 17.3 21.7 24.9 26.3 36.3 ...
 $ MPG.city          : int  25 18 20 19 22 22 19 16 19 16 ...
 $ MPG.highway       : int  31 25 26 26 30 31 28 25 27 25 ...
 $ AirBags           : Factor w/ 3 levels "Driver & Passenger",..: 3 1 2 1 2 2 2 2 2 2 ...
 $ DriveTrain        : Factor w/ 3 levels "4WD","Front",..: 2 2 2 2 3 2 2 3 2 2 ...
 $ Cylinders         : Factor w/ 6 levels "3","4","5","6",..: 2 4 4 4 2 2 4 4 4 5 ...
 $ EngineSize        : num  1.8 3.2 2.8 2.8 3.5 2.2 3.8 5.7 3.8 4.9 ...
 $ Horsepower        : int  140 200 172 172 208 110 170 180 170 200 ...
 $ RPM               : int  6300 5500 5500 5500 5700 5200 4800 4000 4800 4100 ...
 $ Rev.per.mile      : int  2890 2335 2280 2535 2545 2565 1570 1320 1690 1510 ...
 $ Man.trans.avail   : Factor w/ 2 levels "No","Yes": 2 2 2 2 2 1 1 1 1 1 ...
 $ Fuel.tank.capacity: num  13.2 18 16.9 21.1 21.1 16.4 18 23 18.8 18 ...
 $ Passengers        : int  5 5 5 6 4 6 6 6 5 6 ...
 $ Length            : int  177 195 180 193 186 189 200 216 198 206 ...
 $ Wheelbase         : int  102 115 102 106 109 105 111 116 108 114 ...
 $ Width             : int  68 71 67 70 69 69 74 78 73 73 ...
 $ Turn.circle       : int  37 38 37 37 39 41 42 45 41 43 ...
 $ Rear.seat.room    : num  26.5 30 28 31 27 28 30.5 30.5 26.5 35 ...
 $ Luggage.room      : int  11 15 14 17 13 16 17 21 14 18 ...
 $ Weight            : int  2705 3560 3375 3405 3640 2880 3470 4105 3495 3620 ...
 $ Origin            : Factor w/ 2 levels "USA","non-USA": 2 2 2 2 2 1 1 1 1 1 ...
 $ Make              : Factor w/ 93 levels "Acura Integra",..: 1 2 4 3 5 6 7 9 8 10 .

아래의 Histogram, Kernel density curve plot을 보면 Price 는 왼쪽으로 치우쳐있고 오른쪽으로 꼬리가 긴 분포 형태를 띠고 있음을 알 수 있습니다.

> # histogram, kernel density curve
> hist(Cars93$MPG.highway, freq=FALSE, breaks=30)
> 
> hist(Cars93$Price, freq=FALSE, breaks=20)
> lines(density(Cars93$Price), col="red", lty=3)

중심화 경향, 대표값으로 산술평균, 중앙값, 최빈값을 많이 사용하는데요, 모두 나름의 특징과 한계를 가지고 있습니다. 따라서 어느 하나의 통계량이 다른 통계량보다 우수하다거나 더 좋다고 말할 수는 없으며, 분석 목적, 분석 대상 자료의 특성, 업의 특성 등을 종합적으로 감안하여 조심해서 사용, 해석해야 합니다. 그리고 반드시 퍼짐 정도와 분포형태/대칭정도를 나타내는 통계량과 그래프를 병행해서 분석을 해야 왜곡된 해석을 피할 수 있습니다. 아래 개념을 정확히 이해하지 못할 경우 통계를 가지고 사기치는 지능범에게 당하는 수가 있으니 꼭 알아두어야할 기본 개념이 되겠습니다.

(1) 산술평균 (arithmetic mean, average) : mean()

우리들이 일반적으로 "평균"이라고 말할 때 사용하는 것이 바로 "산술평균(arithmetic mean)"입니다. 평균에는 "산술평균" 말고도 "기하평균", "조화평균", "가중평균" 등 여러 종류가 있습니다.

산술평균은 모집단이 정규분포를 띠고 있을 때 가장 적합한 중심화 경향 통계량이라고 할 수 있습니다.

달리 말하면, 분포 형태가 한쪽으로 치우쳐 있다든지, 이상값(outlier)가 있으면 영향을 크게 받으므로 사용에 주의를 요하는 통계량이라고 할 수 있습니다. 이럴 경우에는 정규분포로 변환을 하거나 이상값(outlier)를 제거한 후에 산술평균을 계산하는 것이 바람직한 조치라고 하겠습니다.

> # mean : mean()

> mean(Cars93$Price)
[1] 19.50968

(2) 중앙값 (median) : median()

중앙값은 이상치(outlier)에 덜 민감(robust)하므로 이번 예제처럼 오른쪽으로 긴 꼬리 부분에 초고가(extremely high price)의 차량이 소수 있는 경우의 분포에는 산술평균보다 더 적합한 통계량이라고 할 수 있겠습니다.

> # median : median()
> median(Cars93$Price)
[1] 17.7

(3) 최빈값 (mode) : which.max(table())

최빈값은 연속형 데이터를 가지고 바로 적용해서는 안되며, 사전에 범주형 데이터로 변환을 한 후에, 도수분포표(frequency distribution table)을 작성해서, 도수가 가장 많은 구간(class)을 선정하면 되겠습니다.

> # mode : which.max(table())
> Cars93 <- within(Cars93, {
+   Price_cd = character()
+   Price_cd[Price < 10] = "5_10"
+   Price_cd[Price >= 10 & Price < 15] = "10_15"
+   Price_cd[Price >= 15 & Price < 20] = "15_20"
+   Price_cd[Price >= 20 & Price < 25] = "20_25"
+   Price_cd[Price >= 25 & Price < 30] = "25_30"
+   Price_cd[Price >= 30 & Price < 35] = "30_35"
+   Price_cd[Price >= 35 & Price < 40] = "35_40"
+   Price_cd[Price >= 40 & Price < 45] = "40_45"
+   Price_cd[Price >= 45 & Price < 50] = "45_50"
+   Price_cd[Price >= 50 & Price < 55] = "50_55"
+   Price_cd[Price >= 55 & Price < 60] = "55_60"
+   Price_cd[Price >= 60 ] = "60_65"
+   Price_cd = factor(Price_cd, level=c("5_10", "10_15", "15_20", "20_25", "25_30", "30_35", 
+                                       "35_40", "40_45", "45_50", "50_55", "55)60", "60_65"))
+ })
> 
> table(Cars93$Price_cd)

 5_10 10_15 15_20 20_25 25_30 30_35 35_40 40_45 45_50 50_55 55)60 60_65 
   10    23    28    11     8     6     4     1     1     0     0     1 
> 
> which.max(table(Cars93$Price_cd))
15_20 
    3

3번째 구간인 '15~20' 구간에서 도수가 28개로 가장 많이 나왔으므로, 최빈값은 '15~20' 구간이 되겠습니다.

===================================================== d^_^b =====================================================

아래의 기하평균, CAGR, 조화평균, 가중평균은 우리들이 일상생활에서 산술평균 대비 많이 사용하지는 않습니다만, 산술평균과 많이 헷갈리게 사용하기도 하고, 혹은 산술평균을 적용하면 안되는 상황임에도 아래의 다른 평균을 몰라서 산술평균을 잘못 적용하기도 합니다. 어떠한 상황에서 무슨 대표 통계량을 이용하는지 파악해두면 유용하겠지요?!

(4) 기하평균 (geometric mean) : prod(x)^(1/n), where n = length(x)

기하평균은 인구성장률, 투자이율과 같이 성장률 평균을 산출할 때 사용합니다. 성장률 평균 산출 시 산술평균을 사용하면 안됩니다. 복리 개념의 성장률은 면적의 개념으로 접근을 해야 하므로 기하평균을 사용하게 됩니다.

문제) 작년 이율이 1%, 올해 이율이 5%인 복리정기예금의 2년간 평균 이율은?

~~산술평균 = (1.01 + 1.05)/2 = 1.03, 즉 3.0% (삐이~ 잘못된 계산임)~~
기하평균 = (1.01*1.05)^(1/2) = 1.029806, 즉 2.98% (제대로된 계산임)

위의 문제를 R로 풀어보면 아래와 같습니다.

> # geometric mean
> x <- c(1.01, 1.05) # interest rate of 1st and 2nd year
> prod(x) # 1st year rate x 2nd year rate
[1] 1.0605
> n <- length(x) # length of x vector
> 
> prod(x)^(1/n) # geometric mean
[1] 1.029806

(5) 연평균성장률 (CAGR : componded average growth rate) : (FV/IV)^(1/n)-1

where, IV : initial value of an investment, FV : final value of an investment, n : investment periods

아래 A, B라는 두 회사의 6년간의 매출액과 전년도 대비 성장률을 가지고 예를 들어보겠습니다.

문제) A 회사는 2010년도 부터 해서 100억, 150억, 190억, 250억, 290억, 350억원의 매출을 올렸습니다.

그러면 이 회사의 6년에 걸친 매출액의 연평균성장률은 얼마일까요? 매년 균등하게 몇 %씩 성장했을까요?

~~산술평균(arithmetic mean)으로 계산하면 29.0% 인데요, 이게 맞는 연평균성장률일까요? 땡~ 틀렸습니다.~~

정답은 기하평균 개념을 적용한 CAGR 28.5% 가 되겠습니다.

Company B가 매년 똑같이 28.5%씩 성장한 회사인데요, 2010년도에 100억원에서 시작해서 매년 28.5%씩 성장했더니 2015년에 350억원의 매출이 되어있는 것을 확인할 수 있습니다.

위 예제의 CAGR을 R로 계산해보겠습니다.

> # CAGR (Compound Average Growth Rate)
> IV <- c(100) # initial value of revenue
> FV <- c(350) # final value of reveune
> n <- c(5) # number of year
> CAGR_rev <- (FV/IV)^(1/n) - 1 # CAGR
> CAGR_rev
[1] 0.2847352

(6) 조화평균 (harmonic mean) : 1/mean(1/x)

조화평균(harmonic mean)은 생산성, 효율 등의 평균 산출 시에 사용합니다.

문제) 집에서 학교까지 편도 30km 거리를 갈 때는 시속 30km/h 인 자전거를 타고 갔고, 올 때는 시속 90km/h 인 자동차를 타고 왔을 때 왕복 평균 시속은?

~~위 문제를 산술평균으로 풀어서 만약 (30 + 90)/2 = 60 km/h 라고 한다면, 땡~! 틀린 답입니다.~~

집과 학교 왕복거리는 총 60km 이고, 이때 걸린 시간은 집에서 학교까지 1시간 (30km거리를 시속 30km/h 로 갔으니깐), 학교에서 집으로 오는데 20분 (30km 거리를 시속 90km/h로 왔으니깐) 걸려서 총 1시간 20분이 걸렸습니다. 즉, 올바른 평균 시속은 60/1.333 = 약 45km/h 가 되겠습니다. 이 문제룰 계산할 때 사용하는 평균이 조화평균이 되겠습니다.

> # harmonic mean : 1/mean(1/x)
> km_per_hour <- c(30, 90)
> 
> arithmetic_mean <- mean(km_per_hour)
> arithmetic_mean # wrong answer
[1] 60
> 
> harmonic_mean <- 1/mean(1/km_per_hour)
> harmonic_mean # correct answer
[1] 45

(7) 가중평균 (weighted average) : weighted.mean()

가중평균(weighted average)은 확률, 가중치를 수반하는 평균을 산출할 때 사용합니다.

문제 1) 홍길동씨가 A, B, C 3개 회사 주식에 각각 700만원, 200만원, 100만원씩 총 1,000만원을 투자하여, 각 회사별로 투자 수익율이 15%, 9%, 5% 나왔다. 그러면 홍길동씨의 주식 투자 평균 수익률은?

~~산술평균으로 계산하면 0.097%이지만 이는 틀린 답입니다.~~ 각 회사별로 주식투자한 금액의 비율(여기서는 weight) 가 서로 다르므로, 수익율에다가 투자금의 비율을 가중평균한 값 12.8%가 답이 되겠습니다.

R을 가지고 가중평균을 구할 때는 가중평균의 정의에 맞추어서 계산을 직접해도 되며 (아래 예제의 첫번째 경우), weighted.mean() 함수를 사용(아래 예제의 두번째 경우)해도 되겠습니다.

> # Q1

> # weighted mean
> weighted_earning_rate_1 <- (0.7*0.15 + 0.2*0.09 + 0.1*0.05)/(0.7 + 0.2 + 0.1)
> weighted_earning_rate_1
[1] 0.128
> 
> investment <- data.frame(weight=c(0.7, 0.2, 0.1), earning_rate=c(0.15, 0.09, 0.05))
> weighted_earning_rate_2 <- weighted.mean(investment$earning_rate, investment$weight)
> weighted_earning_rate_2
[1] 0.128

문제 2) 알코올 도수 9%인 와인 200ml와 알코올 도수 21%인 소주 1000ml를 섞어서 와소 폭탄주를 만들었다. 와소 폭탄주의 평균 알코올 도수는?

와소 폭탄주의 평균 알코올 도수 정답은 알코올 도수와 양을 가지고 가중평균으로 구한 19%가 되겠습니다.

R에서는 가중평균의 정의에 따라서 공식에 대해서서 구할 수도 있고 (아래 예의 첫번째 경우), 아니면 weighted.mean() 함수를 사용(아래 예이 두번째 경우)해서 구해도 됩니다.

 > # Q2
> weighted_alcohol_mean_1 <- (200*0.09 + 1000*0.21)/(200+1000)
> weighted_alcohol_mean_1
[1] 0.19
> 
> alcohol <- data.frame(volume=c(200, 1000), alcohol_rate=c(0.09, 0.21))
> weighted_alcohol_mean_2 <- weighted.mean(alcohol$alcohol_rate, alcohol$volume)
> weighted_alcohol_mean_2
[1] 0.19

다음번 포스팅에서는 퍼짐 정도 (dispersion) 에 대한 통계량에 대해서 소개하겠습니다.

중심화 경향과 퍼짐 정도가 다른 두 데이터셋을 표준화하는 방법은 아래의 포스팅을 참고하시기 바랍니다.

☞ R 데이터 변환 (1) 표준화 : z 표준화 변환, [0-1] 변환

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

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

728x90

저작자표시 비영리 변경금지

'R 분석과 프로그래밍 > R 통계분석' 카테고리의 다른 글

R 일변량 연속형 자료 요약 - (3) 분포형태와 대칭정도(distrubution, symmetry) : 왜도(skewness), 첨도(kurtosis), 분위수-분위수(Quantile-Qunantile) (0)	2015.10.12
R 일변량 연속형 자료 요약 - (2) 퍼짐 정도(dispersion) : 분산(variance), 표준편차(sd), 변이계수, 범위(range), IQR, 백분위수(percentile) (0)	2015.10.11
R 이변량 범주형 자료의 정리와 요약 - 분할표(contingency table) : table(), xtabs(), prop.table(), margin.table(), addmargins(), CrossTable() (0)	2015.10.10
R 일변량 범주형 자료의 정리와 요약 - 도수분포표(frequency distribution table) : table(), xtabs() (2)	2015.10.10
R 단일 모집단 분포의 정규성 검정 : shapiro.test(), qqnorm(), qqline() (0)	2015.10.09

Posted by Rfriend

R 이변량 범주형 자료의 정리와 요약 - 분할표(contingency table) : table(), xtabs(), prop.table(), margin.table(), addmargins(), CrossTable()

R 분석과 프로그래밍/R 통계분석 2015. 10. 10. 14:54

자료의 정리와 요약으로

- 일변량 범주형 자료의

도수분포표 (frequency distribution table)

- 이변량 범주형 자료의

분할표 (contingency table)

카이제곱 검정 (chisquare test)

중에서 이번 포스팅에서는 이변량 범주형 자료의 분할표, 카이제곱 검정에 대해서 알아보겠습니다.

참고로, 자료의 형태(질적/범주형변수 vs. 양적변수)에 따라서, 또 양적변수는 연속형이냐 이산형이냐에 따라서 자료를 정리하고 요약하는 방법, 확률모형 설정, 추정과 검정, 모형 개발 등이 달라지는데요, 아래 표를 참고하시기 바랍니다.

[ 자료의 형태 구분 ]

형태 (Type)	순서 여부 (ordinal yes/no) 연속성 여부 (continuous yes/no)	예시 (example)
질적 변수 (Qualitative variable) or 범주형 변수 (Categorical variable)	no : 명목형 (nominal)	이름, 성별, 주소, 학력, 전화번호 등
	yes : 순서형 (ordinal)	학년
양적 변수 (Quantitative variable)	yes : 연속형 변수 (continuous variable)	키, 몸무게, 온도, 습도 등
양적 변수 (Quantitative variable)	no : 이산형 변수 (descrete variable)	나이

분석 기법은 변수가 몇개냐, 변수의 형태가 무엇이냐(범주형 변수와 연속형 변수의 2개로 구분)해서 간략하게 표로 제시해보면 아래와 같습니다.

일변량 (one variable)

자료 형태

분석 기법

범주형 변수

(categorical variable)

<표(table)를 이용한 정리와 요약 >

도수분포표 (frequency distribution table)

연속형 변수

(continuous variable)

< 통계량(statistics)을 이용한 요약 >

- 중심 경향 : 평균 (mean), 중앙값 (median), 최빈값 (mode)

- 퍼짐 정도 : 분산 (variance), 표준편차(standard deviation),

범위 (range), 사분위수(IQR: Inter Quantile Rnage)

- 치우침 정도 : 왜도 (skewness), 첨도 (kurtosis)

※ 그래프 분석 병행 매우 중요

이변량 (two variables)

Y 변수

X 변수

범주형 변수

(categorical variable)

연속형 변수

(continuous variable)

범주형 변수

(categorical variable)

분할표 (contingency table)

카이제곱 검정 (chi-square test)

t-Test

ANOVA

연속형 변수

(continuous variable)

의사결정나무 (decision tree)

로지스틱회귀분석 (rogistic regression)

상관분석 (correlation analysis)

회귀분석 (regression analysis)

※ 그래프 분석 병행 매우 중요

일변량 자료 분석이 변수 자체의 분포, 형태에 대해 관심을 가진다면, 이변량 자료 분석의 경우는 변수와 변수간의 관계에 특히 관심을 가지게 됩니다.

R 예제에 사용한 데이터는 MASS 패키지에 내장된 Cars93 데이터프레임의 차종형태(Type), 실린더(Cylinders), 트렁크공간크기(Luggage.room)의 변수들입니다.

Type, Cylinders 변수는 categorical data 이어서 교차표(crosstabulation table) 분석할 때 그대로 사용하면 되며, Luggage.room 변수는 discrete data (integer) 이므로 categorical data로 변환해서 분할표 분석할 때 사용하겠습니다.

> library(MASS)
> str(Cars93)
'data.frame':	93 obs. of  27 variables:
 $ Manufacturer      : Factor w/ 32 levels "Acura","Audi",..: 1 1 2 2 3 4 4 4 4 5 ...
 $ Model             : Factor w/ 93 levels "100","190E","240",..: 49 56 9 1 6 24 54 74 73 35 ...
 $ Type              : Factor w/ 6 levels "Compact","Large",..: 4 3 1 3 3 3 2 2 3 2 ...
 $ Min.Price         : num  12.9 29.2 25.9 30.8 23.7 14.2 19.9 22.6 26.3 33 ...
 $ Price             : num  15.9 33.9 29.1 37.7 30 15.7 20.8 23.7 26.3 34.7 ...
 $ Max.Price         : num  18.8 38.7 32.3 44.6 36.2 17.3 21.7 24.9 26.3 36.3 ...
 $ MPG.city          : int  25 18 20 19 22 22 19 16 19 16 ...
 $ MPG.highway       : int  31 25 26 26 30 31 28 25 27 25 ...
 $ AirBags           : Factor w/ 3 levels "Driver & Passenger",..: 3 1 2 1 2 2 2 2 2 2 ...
 $ DriveTrain        : Factor w/ 3 levels "4WD","Front",..: 2 2 2 2 3 2 2 3 2 2 ...
 $ Cylinders         : Factor w/ 6 levels "3","4","5","6",..: 2 4 4 4 2 2 4 4 4 5 ...
 $ EngineSize        : num  1.8 3.2 2.8 2.8 3.5 2.2 3.8 5.7 3.8 4.9 ...
 $ Horsepower        : int  140 200 172 172 208 110 170 180 170 200 ...
 $ RPM               : int  6300 5500 5500 5500 5700 5200 4800 4000 4800 4100 ...
 $ Rev.per.mile      : int  2890 2335 2280 2535 2545 2565 1570 1320 1690 1510 ...
 $ Man.trans.avail   : Factor w/ 2 levels "No","Yes": 2 2 2 2 2 1 1 1 1 1 ...
 $ Fuel.tank.capacity: num  13.2 18 16.9 21.1 21.1 16.4 18 23 18.8 18 ...
 $ Passengers        : int  5 5 5 6 4 6 6 6 5 6 ...
 $ Length            : int  177 195 180 193 186 189 200 216 198 206 ...
 $ Wheelbase         : int  102 115 102 106 109 105 111 116 108 114 ...
 $ Width             : int  68 71 67 70 69 69 74 78 73 73 ...
 $ Turn.circle       : int  37 38 37 37 39 41 42 45 41 43 ...
 $ Rear.seat.room    : num  26.5 30 28 31 27 28 30.5 30.5 26.5 35 ...
 $ Luggage.room      : int  11 15 14 17 13 16 17 21 14 18 ...
 $ Weight            : int  2705 3560 3375 3405 3640 2880 3470 4105 3495 3620 ...
 $ Origin            : Factor w/ 2 levels "USA","non-USA": 2 2 2 2 2 1 1 1 1 1 ...
 $ Make              : Factor w/ 93 levels "Acura Integra",..: 1 2 4 3 5 6 7 9 8 10 ..

R로 실습을 해보겠습니다.

(1) 이변량 분할표 (contingency table) : table(), xtabs()

> ##----------------------------------------------------------
> ## 이변량 분할표(contingency table) : table(), xtabs()
> ##----------------------------------------------------------
>

> # contingency table : table(x, y)
> Car_table_3 <- with(Cars93, table(Type, Cylinders))
> Car_table_3
         Cylinders
Type       3  4  5  6  8 rotary
  Compact  0 15  0  1  0      0
  Large    0  0  0  7  4      0
  Midsize  0  7  1 12  2      0
  Small    3 18  0  0  0      0
  Sporty   0  8  0  4  1      1
  Van      0  1  1  7  0      0
> 
> 
> # contingency table : xtabs(~ x + y, data)
> Car_table_4 <- xtabs(~ Type + Cylinders, data=Cars93)
> Car_table_4
         Cylinders
Type       3  4  5  6  8 rotary
  Compact  0 15  0  1  0      0
  Large    0  0  0  7  4      0
  Midsize  0  7  1 12  2      0
  Small    3 18  0  0  0      0
  Sporty   0  8  0  4  1      1
  Van      0  1  1  7  0      0

참고로, table 형태의 데이터를 가지고 xtabs() 함수를 사용해서 cross tabulation 계산할 때의 예를 아래에 들어보겠습니다. 데이터 형태가 위의 Cars93과는 다르다보니 xtabs() 함수 쓸 때 '~' 앞에 도수가 들어있는 변수를 지정해주는 것이 다릅니다. 아래에 UC버클리 대학교 입학 관련 데이터에의 예의 제일 아래에 보면 xtabs(Freq ~ Gender + Admit, data=UCBAdmissions.df) 라고 했지요? 위의 Cars93의 예의 경우는 xtabs( ~ Type + Cylinders, data=Cars93)이라고 했구요. 데이터가 어떻게 생겼느냐에 따라서 알맞는 arguments를 사용하면 되겠습니다.

> # getting data of table format > data(UCBAdmissions) > str(UCBAdmissions) table [1:2, 1:2, 1:6] 512 313 89 19 353 207 17 8 120 205 ... - attr(*, "dimnames")=List of 3 ..$ Admit : chr [1:2] "Admitted" "Rejected" ..$ Gender: chr [1:2] "Male" "Female" ..$ Dept : chr [1:6] "A" "B" "C" "D" ... > UCBAdmissions , , Dept = A Gender Admit Male Female Admitted 512 89 Rejected 313 19 , , Dept = B Gender Admit Male Female Admitted 353 17 Rejected 207 8 , , Dept = C Gender Admit Male Female Admitted 120 202 Rejected 205 391 , , Dept = D Gender Admit Male Female Admitted 138 131 Rejected 279 244 , , Dept = E Gender Admit Male Female Admitted 53 94 Rejected 138 299 , , Dept = F Gender Admit Male Female Admitted 22 24 Rejected 351 317 > > # transforming to dataframe > UCBAdmissions.df <- as.data.frame(UCBAdmissions) > UCBAdmissions.df Admit Gender Dept Freq 1 Admitted Male A 512 2 Rejected Male A 313 3 Admitted Female A 89 4 Rejected Female A 19 5 Admitted Male B 353 6 Rejected Male B 207 7 Admitted Female B 17 8 Rejected Female B 8 9 Admitted Male C 120 10 Rejected Male C 205 11 Admitted Female C 202 12 Rejected Female C 391 13 Admitted Male D 138 14 Rejected Male D 279 15 Admitted Female D 131 16 Rejected Female D 244 17 Admitted Male E 53 18 Rejected Male E 138 19 Admitted Female E 94 20 Rejected Female E 299 21 Admitted Male F 22 22 Rejected Male F 351 23 Admitted Female F 24 24 Rejected Female F 317 > > # cross tabulation of transformed dataset > xtabs(Freq ~ Gender + Admit, data=UCBAdmissions.df) Admit Gender Admitted Rejected Male 1198 1493 Female 557 1278

(2) 상대 분할표 (relative frequency contingency table) : prop.table()

> # relative frequency distribution table : prop.table()
> options("digit" = 3) # decimal point setting
> Prop_Car_table_3 <- prop.table(Car_table_3)
> Prop_Car_table_3
         Cylinders
Type               3          4          5          6          8     rotary
  Compact 0.00000000 0.16129032 0.00000000 0.01075269 0.00000000 0.00000000
  Large   0.00000000 0.00000000 0.00000000 0.07526882 0.04301075 0.00000000
  Midsize 0.00000000 0.07526882 0.01075269 0.12903226 0.02150538 0.00000000
  Small   0.03225806 0.19354839 0.00000000 0.00000000 0.00000000 0.00000000
  Sporty  0.00000000 0.08602151 0.00000000 0.04301075 0.01075269 0.01075269
  Van     0.00000000 0.01075269 0.01075269 0.07526882 0.00000000 0.00000000

(3) 한계분포 (marginal distribution) : margin.table(table, margin=1 or 2)

> margin.table(Car_table_3, margin=1) # margin=1 for row sum
Type
Compact   Large Midsize   Small  Sporty     Van 
     16      11      22      21      14       9 
> margin.table(Prop_Car_table_3, margin=1) # margin=1 for row sum
Type
   Compact      Large    Midsize      Small     Sporty        Van 
0.17204301 0.11827957 0.23655914 0.22580645 0.15053763 0.09677419 
> 
> margin.table(Car_table_3, margin=2) # margin=2 for column sum
Cylinders
     3      4      5      6      8 rotary 
     3     49      2     31      7      1 
> margin.table(Prop_Car_table_3, margin=2) # margin=2 for column sum
Cylinders
         3          4          5          6          8     rotary 
0.03225806 0.52688172 0.02150538 0.33333333 0.07526882 0.01075269

(4) 분할표에 한계분포 추가 (add marginal distribution at the table) addmargins()

> # marginal distribution at the table : addmargins(table, margin=1 or 2)
> addmargins(Car_table_3) # row sum and column sum at the table
         Cylinders
Type       3  4  5  6  8 rotary Sum
  Compact  0 15  0  1  0      0  16
  Large    0  0  0  7  4      0  11
  Midsize  0  7  1 12  2      0  22
  Small    3 18  0  0  0      0  21
  Sporty   0  8  0  4  1      1  14
  Van      0  1  1  7  0      0   9
  Sum      3 49  2 31  7      1  93
> addmargins(Prop_Car_table_3) # row sum and column sum at the table
         Cylinders
Type               3          4          5          6          8     rotary        Sum
  Compact 0.00000000 0.16129032 0.00000000 0.01075269 0.00000000 0.00000000 0.17204301
  Large   0.00000000 0.00000000 0.00000000 0.07526882 0.04301075 0.00000000 0.11827957
  Midsize 0.00000000 0.07526882 0.01075269 0.12903226 0.02150538 0.00000000 0.23655914
  Small   0.03225806 0.19354839 0.00000000 0.00000000 0.00000000 0.00000000 0.22580645
  Sporty  0.00000000 0.08602151 0.00000000 0.04301075 0.01075269 0.01075269 0.15053763
  Van     0.00000000 0.01075269 0.01075269 0.07526882 0.00000000 0.00000000 0.09677419
  Sum     0.03225806 0.52688172 0.02150538 0.33333333 0.07526882 0.01075269 1.00000000
> 
>

> addmargins(Car_table_3, margin=1) # margin=1 for row sum at the table
         Cylinders
Type       3  4  5  6  8 rotary
  Compact  0 15  0  1  0      0
  Large    0  0  0  7  4      0
  Midsize  0  7  1 12  2      0
  Small    3 18  0  0  0      0
  Sporty   0  8  0  4  1      1
  Van      0  1  1  7  0      0
  Sum      3 49  2 31  7      1
> addmargins(Prop_Car_table_3, margin=1) # margin=1 for row sum at the table
         Cylinders
Type               3          4          5          6          8     rotary
  Compact 0.00000000 0.16129032 0.00000000 0.01075269 0.00000000 0.00000000
  Large   0.00000000 0.00000000 0.00000000 0.07526882 0.04301075 0.00000000
  Midsize 0.00000000 0.07526882 0.01075269 0.12903226 0.02150538 0.00000000
  Small   0.03225806 0.19354839 0.00000000 0.00000000 0.00000000 0.00000000
  Sporty  0.00000000 0.08602151 0.00000000 0.04301075 0.01075269 0.01075269
  Van     0.00000000 0.01075269 0.01075269 0.07526882 0.00000000 0.00000000
  Sum     0.03225806 0.52688172 0.02150538 0.33333333 0.07526882 0.01075269
>  
>

> addmargins(Car_table_3, margin=2) # margin=2 for column sum at the table
         Cylinders
Type       3  4  5  6  8 rotary Sum
  Compact  0 15  0  1  0      0  16
  Large    0  0  0  7  4      0  11
  Midsize  0  7  1 12  2      0  22
  Small    3 18  0  0  0      0  21
  Sporty   0  8  0  4  1      1  14
  Van      0  1  1  7  0      0   9
> addmargins(Prop_Car_table_3, margin=2) # margin=2 for column sum at the table
         Cylinders
Type               3          4          5          6          8     rotary        Sum
  Compact 0.00000000 0.16129032 0.00000000 0.01075269 0.00000000 0.00000000 0.17204301
  Large   0.00000000 0.00000000 0.00000000 0.07526882 0.04301075 0.00000000 0.11827957
  Midsize 0.00000000 0.07526882 0.01075269 0.12903226 0.02150538 0.00000000 0.23655914
  Small   0.03225806 0.19354839 0.00000000 0.00000000 0.00000000 0.00000000 0.22580645
  Sporty  0.00000000 0.08602151 0.00000000 0.04301075 0.01075269 0.01075269 0.15053763
  Van     0.00000000 0.01075269 0.01075269 0.07526882 0.00000000 0.00000000 0.09677419

margin.table()은 row 혹은 column 별로 sum을 따로 구한데 반해, addmargins() 는 contingency table은 고정시켜놓고 거기에 row 혹은 column sum을 추가한 것이 다릅니다.

(5) NA를 분할표에 포함시키고자 할 경우 (including NA to contingency table) : useNA = "ifany"

> ##------
> # useNA = "ifany" option
> # count of missing value
> dim(Cars93)
[1] 93 27
> sum(is.na(Cars93$Type))
[1] 0
> sum(is.na(Cars93$Cylinders))
[1] 0
> sum(is.na(Cars93$Luggage.room)) # number of NA = 11
[1] 11
> summary(Cars93$Luggage.room)
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max.    NA's 
   6.00   12.00   14.00   13.89   15.00   22.00      11 
> 
> # useNA = "ifany" option
> addmargins(with(Cars93, table((Luggage.room_big = Luggage.room >= 14), Type)))
       Type
        Compact Large Midsize Small Sporty Van Sum
  FALSE       4     0       2    18     10   0  34
  TRUE       12    11      20     3      2   0  48
  Sum        16    11      22    21     12   0  82
> addmargins(with(Cars93, table((Luggage.room_big = Luggage.room >= 14), Type, useNA="ifany")))
       Type
        Compact Large Midsize Small Sporty Van Sum
  FALSE       4     0       2    18     10   0  34
  TRUE       12    11      20     3      2   0  48
  <NA>        0     0       0     0      2   9  11
  Sum        16    11      22    21     14   9  93

(6) 카이제곱 검정 (chi-square test) : gmodels package, CrossTable(x,y, expected=T, chisq=T)

> install.packages("gmodels")
Installing package into ‘C:/Users/user/Documents/R/win-library/3.2’
(as ‘lib’ is unspecified)
trying URL 'http://cran.rstudio.com/bin/windows/contrib/3.2/gmodels_2.16.2.zip'
Content type 'application/zip' length 73944 bytes (72 KB)
downloaded 72 KB

package ‘gmodels’ successfully unpacked and MD5 sums checked

The downloaded binary packages are in
	C:\Users\user\AppData\Local\Temp\RtmpSUTL1V\downloaded_packages
> library(gmodels)

> with(Cars93, CrossTable(Type, Cylinders, expected=TRUE, chisq=TRUE))

 
   Cell Contents
|-------------------------|
|                       N |
|              Expected N |
| Chi-square contribution |
|           N / Row Total |
|           N / Col Total |
|         N / Table Total |
|-------------------------|

 
Total Observations in Table:  93 

 
             | Cylinders 
        Type |         3 |         4 |         5 |         6 |         8 |    rotary | Row Total | 
-------------|-----------|-----------|-----------|-----------|-----------|-----------|-----------|
     Compact |         0 |        15 |         0 |         1 |         0 |         0 |        16 | 
             |     0.516 |     8.430 |     0.344 |     5.333 |     1.204 |     0.172 |           | 
             |     0.516 |     5.120 |     0.344 |     3.521 |     1.204 |     0.172 |           | 
             |     0.000 |     0.938 |     0.000 |     0.062 |     0.000 |     0.000 |     0.172 | 
             |     0.000 |     0.306 |     0.000 |     0.032 |     0.000 |     0.000 |           | 
             |     0.000 |     0.161 |     0.000 |     0.011 |     0.000 |     0.000 |           | 
-------------|-----------|-----------|-----------|-----------|-----------|-----------|-----------|
       Large |         0 |         0 |         0 |         7 |         4 |         0 |        11 | 
             |     0.355 |     5.796 |     0.237 |     3.667 |     0.828 |     0.118 |           | 
             |     0.355 |     5.796 |     0.237 |     3.030 |    12.153 |     0.118 |           | 
             |     0.000 |     0.000 |     0.000 |     0.636 |     0.364 |     0.000 |     0.118 | 
             |     0.000 |     0.000 |     0.000 |     0.226 |     0.571 |     0.000 |           | 
             |     0.000 |     0.000 |     0.000 |     0.075 |     0.043 |     0.000 |           | 
-------------|-----------|-----------|-----------|-----------|-----------|-----------|-----------|
     Midsize |         0 |         7 |         1 |        12 |         2 |         0 |        22 | 
             |     0.710 |    11.591 |     0.473 |     7.333 |     1.656 |     0.237 |           | 
             |     0.710 |     1.819 |     0.587 |     2.970 |     0.071 |     0.237 |           | 
             |     0.000 |     0.318 |     0.045 |     0.545 |     0.091 |     0.000 |     0.237 | 
             |     0.000 |     0.143 |     0.500 |     0.387 |     0.286 |     0.000 |           | 
             |     0.000 |     0.075 |     0.011 |     0.129 |     0.022 |     0.000 |           | 
-------------|-----------|-----------|-----------|-----------|-----------|-----------|-----------|
       Small |         3 |        18 |         0 |         0 |         0 |         0 |        21 | 
             |     0.677 |    11.065 |     0.452 |     7.000 |     1.581 |     0.226 |           | 
             |     7.963 |     4.347 |     0.452 |     7.000 |     1.581 |     0.226 |           | 
             |     0.143 |     0.857 |     0.000 |     0.000 |     0.000 |     0.000 |     0.226 | 
             |     1.000 |     0.367 |     0.000 |     0.000 |     0.000 |     0.000 |           | 
             |     0.032 |     0.194 |     0.000 |     0.000 |     0.000 |     0.000 |           | 
-------------|-----------|-----------|-----------|-----------|-----------|-----------|-----------|
      Sporty |         0 |         8 |         0 |         4 |         1 |         1 |        14 | 
             |     0.452 |     7.376 |     0.301 |     4.667 |     1.054 |     0.151 |           | 
             |     0.452 |     0.053 |     0.301 |     0.095 |     0.003 |     4.793 |           | 
             |     0.000 |     0.571 |     0.000 |     0.286 |     0.071 |     0.071 |     0.151 | 
             |     0.000 |     0.163 |     0.000 |     0.129 |     0.143 |     1.000 |           | 
             |     0.000 |     0.086 |     0.000 |     0.043 |     0.011 |     0.011 |           | 
-------------|-----------|-----------|-----------|-----------|-----------|-----------|-----------|
         Van |         0 |         1 |         1 |         7 |         0 |         0 |         9 | 
             |     0.290 |     4.742 |     0.194 |     3.000 |     0.677 |     0.097 |           | 
             |     0.290 |     2.953 |     3.360 |     5.333 |     0.677 |     0.097 |           | 
             |     0.000 |     0.111 |     0.111 |     0.778 |     0.000 |     0.000 |     0.097 | 
             |     0.000 |     0.020 |     0.500 |     0.226 |     0.000 |     0.000 |           | 
             |     0.000 |     0.011 |     0.011 |     0.075 |     0.000 |     0.000 |           | 
-------------|-----------|-----------|-----------|-----------|-----------|-----------|-----------|
Column Total |         3 |        49 |         2 |        31 |         7 |         1 |        93 | 
             |     0.032 |     0.527 |     0.022 |     0.333 |     0.075 |     0.011 |           | 
-------------|-----------|-----------|-----------|-----------|-----------|-----------|-----------|

 
Statistics for All Table Factors


Pearson's Chi-squared test 
------------------------------------------------------------
Chi^2 =  78.93491     d.f. =  25     p =  1.674244e-07

P-value 가 매우 작으므로 차종(Type)과 Cylinder개수(Cylinders)는 5% 유의수준(significnace level)에서 서로 독립이 아니라(귀무가설 기각, 대립가설 채택)고 판단할 수 있겠습니다.

(7) Mosaic plot for categorical data : vcd package, mosaic()


> # mosaic plot : vcd package, mosaic()

> library(vcd)
필요한 패키지를 로딩중입니다: grid
> Car_table_3 <- with(Cars93, table(Type, Cylinders))
> Car_table_3
         Cylinders
Type       3  4  5  6  8 rotary
  Compact  0 15  0  1  0      0
  Large    0  0  0  7  4      0
  Midsize  0  7  1 12  2      0
  Small    3 18  0  0  0      0
  Sporty   0  8  0  4  1      1
  Van      0  1  1  7  0      0
> mosaic(Car_table_3, 
+        gp=gpar(fill=c("red", "blue")), 
+        direction="v", 
+        main="Mosaic plot of Car Type & Cylinders")

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

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

728x90

저작자표시 비영리 변경금지

'R 분석과 프로그래밍 > R 통계분석' 카테고리의 다른 글

R 일변량 연속형 자료 요약 - (2) 퍼짐 정도(dispersion) : 분산(variance), 표준편차(sd), 변이계수, 범위(range), IQR, 백분위수(percentile) (0)	2015.10.11
R 일변량 연속형 자료 요약 - (1) 중심화 경향 : 산술평균(mean), 중앙값(median), 최빈값, 기하평균, CAGR, 조화평균, 가중평균 (2)	2015.10.11
R 일변량 범주형 자료의 정리와 요약 - 도수분포표(frequency distribution table) : table(), xtabs() (2)	2015.10.10
R 단일 모집단 분포의 정규성 검정 : shapiro.test(), qqnorm(), qqline() (0)	2015.10.09
R 단일 모집단 중심에 대한 비모수 검정 : wilcox.test() (0)	2015.10.04

Posted by Rfriend

R 일변량 범주형 자료의 정리와 요약 - 도수분포표(frequency distribution table) : table(), xtabs()

R 분석과 프로그래밍/R 통계분석 2015. 10. 10. 01:25

자료의 정리와 요약으로

- 일변량 범주형 자료의

도수분포표 (frequency distribution table)

- 이변량 범주형 자료의

분할표 (contingency table)

카이제곱 검정 (chisquare test)

중에서 이번 포스팅에서는 일변량 범주형 자료에 대패 도수분포표(frequency distribution table), 상대도수분포표(relative frequency distribution table)를 사용하여 원자료(raw data)의 정리와 요약에 대해서 알아보도록 하겠습니다.

원자료(raw data)를 수집하고 나면 많은 경우 자료의 양이 너무 방대하기 때문에 모집단의 특성, 형태를 파악하기 위해 자료를 정리, 요약할 필요가 생깁니다.

영화 Matrix에서 데이터가 어마무시하게 쏟아지는 가운데서 뭔가 의미있는 정보 (agency 출현?)뽑아내는 Operator 'TANK' 기억나는지요? 눈이 빠져라 데이터 홍수 속에서 신호와 소음을 구분하느라 용쓰는게 안쓰러워보이지요? ^^;

* 이미지 출처: http://matrixcommunity.org

자료의 형태(질적/범주형변수 vs. 양적변수)에 따라서, 또 양적변수는 연속형이냐 이산형이냐에 따라서 자료를 정리하고 요약하는 방법, 확률모형 설정, 추정과 검정, 모형 개발 등이 달라집니다.

[ 자료의 형태 구분 ]

형태 (Type)	순서 여부 (ordinal yes/no) 연속성 여부 (continuous yes/no)	예시 (example)
질적 변수 (Qualitative variable) or 범주형 변수 (Categorical variable)	no : 명목형 (nominal)	이름, 성별, 주소, 학력, 전화번호 등
	yes : 순서형 (ordinal)	학년
양적 변수 (Quantitative variable)	yes : 연속형 변수 (continuous variable)	키, 몸무게, 온도, 습도 등
양적 변수 (Quantitative variable)	no : 이산형 변수 (descrete variable)	나이

분석 기법은 변수가 몇개냐, 변수의 형태가 무엇이냐(범주형 변수와 연속형 변수의 2개로 구분)해서 간략하게 표로 제시해보면 아래와 같습니다.

일변량 (one variable)

자료 형태

분석 기법

범주형 변수

(categorical variable)

<표(table)를 이용한 정리와 요약 >

도수분포표 (frequency distribution table)

연속형 변수

(continuous variable)

< 통계량(statistics)을 이용한 요약 >

- 중심 경향 : 평균 (mean), 중앙값 (median), 최빈값 (mode)

- 퍼짐 정도 : 분산 (variance), 표준편차(standard deviation),

범위 (range), 사분위수(IQR: Inter Quantile Rnage)

- 치우침 정도 : 왜도 (skewness), 첨도 (kurtosis)

※ 그래프 분석 병행 매우 중요

이변량 (two variables)

Y 변수

X 변수

범주형 변수

(categorical variable)

연속형 변수

(continuous variable)

범주형 변수

(categorical variable)

분할표 (contingency table)

카이제곱 검정 (chi-square test)

t-Test

ANOVA

연속형 변수

(continuous variable)

의사결정나무 (decision tree)

로지스틱회귀분석 (rogistic regression)

상관분석 (correlation analysis)

회귀분석 (regression analysis)

※ 그래프 분석 병행 매우 중요

그럼, R의 table(), xtabs() 함수를 이용하여 도수분포표, 상대도수분포표로 일변량(one variable) 범주형 자료(categorical variable)에 대한 정리와 요약을 해보겠습니다.

MASS 패키지에 내장된 Cars93 데이터프레임의 차종형태(Type) 변수를 분석해 보겠습니다.

> library(MASS)
> str(Cars93)
'data.frame':	93 obs. of  27 variables:
 $ Manufacturer      : Factor w/ 32 levels "Acura","Audi",..: 1 1 2 2 3 4 4 4 4 5 ...
 $ Model             : Factor w/ 93 levels "100","190E","240",..: 49 56 9 1 6 24 54 74 73 35 ...
 $ Type              : Factor w/ 6 levels "Compact","Large",..: 4 3 1 3 3 3 2 2 3 2 ...
 $ Min.Price         : num  12.9 29.2 25.9 30.8 23.7 14.2 19.9 22.6 26.3 33 ...
 $ Price             : num  15.9 33.9 29.1 37.7 30 15.7 20.8 23.7 26.3 34.7 ...
 $ Max.Price         : num  18.8 38.7 32.3 44.6 36.2 17.3 21.7 24.9 26.3 36.3 ...
 $ MPG.city          : int  25 18 20 19 22 22 19 16 19 16 ...
 $ MPG.highway       : int  31 25 26 26 30 31 28 25 27 25 ...
 $ AirBags           : Factor w/ 3 levels "Driver & Passenger",..: 3 1 2 1 2 2 2 2 2 2 ...
 $ DriveTrain        : Factor w/ 3 levels "4WD","Front",..: 2 2 2 2 3 2 2 3 2 2 ...
 $ Cylinders         : Factor w/ 6 levels "3","4","5","6",..: 2 4 4 4 2 2 4 4 4 5 ...
 $ EngineSize        : num  1.8 3.2 2.8 2.8 3.5 2.2 3.8 5.7 3.8 4.9 ...
 $ Horsepower        : int  140 200 172 172 208 110 170 180 170 200 ...
 $ RPM               : int  6300 5500 5500 5500 5700 5200 4800 4000 4800 4100 ...
 $ Rev.per.mile      : int  2890 2335 2280 2535 2545 2565 1570 1320 1690 1510 ...
 $ Man.trans.avail   : Factor w/ 2 levels "No","Yes": 2 2 2 2 2 1 1 1 1 1 ...
 $ Fuel.tank.capacity: num  13.2 18 16.9 21.1 21.1 16.4 18 23 18.8 18 ...
 $ Passengers        : int  5 5 5 6 4 6 6 6 5 6 ...
 $ Length            : int  177 195 180 193 186 189 200 216 198 206 ...
 $ Wheelbase         : int  102 115 102 106 109 105 111 116 108 114 ...
 $ Width             : int  68 71 67 70 69 69 74 78 73 73 ...
 $ Turn.circle       : int  37 38 37 37 39 41 42 45 41 43 ...
 $ Rear.seat.room    : num  26.5 30 28 31 27 28 30.5 30.5 26.5 35 ...
 $ Luggage.room      : int  11 15 14 17 13 16 17 21 14 18 ...
 $ Weight            : int  2705 3560 3375 3405 3640 2880 3470 4105 3495 3620 ...
 $ Origin            : Factor w/ 2 levels "USA","non-USA": 2 2 2 2 2 1 1 1 1 1 ...
 $ Make              : Factor w/ 93 levels "Acura Integra",..: 1 2 4 3 5 6 7 9 8 10 ...

(1) 도수분포표(frequency distribution table) : table(), xtabs()

> ##----------------------------------------------------------
> ## 일변량 도수분포표(frequency distribution table) : table(), xtabs()
> ##----------------------------------------------------------
>

> # frequency distribution table : table()
> Car_table_1 <- with(Cars93, table(Type))
> Car_table_1
Type
Compact   Large Midsize   Small  Sporty     Van 
     16      11      22      21      14       9 
> 
> # frequency distribution table : xtabs()
> Car_table_2 <- xtabs(~Type, data=Cars93)
> Car_table_2
Type
Compact   Large Midsize   Small  Sporty     Van 
     16      11      22      21      14       9

table()은 데이터 프레임을 할당할 수 있는 없어서 with() 함수를 사용했습니다. (attach() 함수를 사용해도 됨).

반면에 xtabs()는 xtabs(formula, data) 형식으로서 뒤에 데이터 프레임을 명시적으로 할당하는 문법입니다.

(2) 상대도수분포표(relative frequency distribution table) : prop.table()

> options("digits" = 2) # 소수점 자리수 설정 (decimal point setting)
> prop.table(Car_table_1)
Type
Compact   Large Midsize   Small  Sporty     Van 
  0.172   0.118   0.237   0.226   0.151   0.097

(3) Bar plot : barplot()

> # bar plot
> barplot(Car_table_1, 
+         main = "Bar Plot of Car Type", 
+         xlab = "Car Type", ylab = "conuts", 
+         col = c("yellow"))

아무래도 그래프로 보는 것이 표로 보는 것보다는 직관적인 이해, 가독성이 우수합니다.

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

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

728x90

저작자표시 비영리 변경금지

'R 분석과 프로그래밍 > R 통계분석' 카테고리의 다른 글

R 일변량 연속형 자료 요약 - (1) 중심화 경향 : 산술평균(mean), 중앙값(median), 최빈값, 기하평균, CAGR, 조화평균, 가중평균 (2)	2015.10.11
R 이변량 범주형 자료의 정리와 요약 - 분할표(contingency table) : table(), xtabs(), prop.table(), margin.table(), addmargins(), CrossTable() (0)	2015.10.10
R 단일 모집단 분포의 정규성 검정 : shapiro.test(), qqnorm(), qqline() (0)	2015.10.09
R 단일 모집단 중심에 대한 비모수 검정 : wilcox.test() (0)	2015.10.04
R 단일 모집단의 모비율에 대한 신뢰구간 추정과 검정(population proportion test) : prop.test() (0)	2015.10.03

Posted by Rfriend

R 단일 모집단 분포의 정규성 검정 : shapiro.test(), qqnorm(), qqline()

R 분석과 프로그래밍/R 통계분석 2015. 10. 9. 13:07

통계적 검정 (statistical testing) 은 모집단의 모수 또는 분포 형태에 대한 추정에 대해 그것이 옳은지 그른지를 임의로 추출한 표본으로부터 통계량을 측정하여 판단하는 통계적 절차를 말합니다.

단일 모집단에 대한 통계적 추론 (추정과 검정) 과 관련하여

[표본이 크고 정규성 충족 시]

- 단일 모집단의 모평균에 대한 신뢰구간 추정과 검정

: t.test()

- 단일 모집단의 모분산에 대한 신뢰구간 추정과 검정

: chi-square test

- 단일 모집단의 모비율에 대한 신뢰구간 추정과 검정

: prop.test()

[정규성 미충족 시]

- 단일 모집단 중심에 대한 비모수 검정

(Nonparametric test on one sample with median)

: wilcox.test()

[정규성 여부 검정]

- 단일 모집단 분포의 정규성 검정

: shapiro.test(), qqnorm(), qqline()

을 차례로 살펴보겠습니다.

지난번 포스팅의 정규분포 형태를 띠는 '단일 모집단의 모평균, 모분산, 모비율에 대한 신뢰구간 추정과 검정'과 정규분포 형태를 띠지 않는 모집단에 대한 비모수 검정을 알아보았습니다. 아래의 표처럼 모집단이 정규성을 만족하느냐 아니냐에 따라서 모수 검정과 비모수 검정으로 분석 기법이 달라지게 됩니다.

[ 모수 검정과 비모수 검정 비교표 ]

구분	모수 검정 (Parametric Test)	비모수 검정 (Nonparametric Test)
When to use	정규분포 가정 만족 시 (normal distribution)	정규분포 가정 불충족 시, 혹은 모집단이 어떤 분포를 따르는지 모를 때 (non-normal distribution, or un-know distribution, or very small sample size, or rankded data)
Statisctic	평균 (mean)	중앙값 (median)
1 sample	1 sample t-test	1 sample Wilcoxon signed rank test
2 samples	2 sample t-test	Mann-Whitney test
2 samples	paired 2-sample t-test	Wilcoxon signed rank test
more than 2 samples	one-way ANOVA	Kruskal-Wallis test

회귀분석, 다변량 통계분석을 할 때도 독립변수의 정규성 가정을 충족시켜야 하며, 이상치(outlier)나 영향값(influencer)가 있을 경우 모형이 매우 민감하게 반응해 모집단을 대표하는 일반화된 모형을 만들 수 없게 됩니다. 따라서 정규성 검정 후에 이상치, 영향치가 있다고 판단되면 이를 제거하는 것이 중요하다고 하겠습니다.

그럼, 이번 포스팅에서는 정규분포 형태를 띠는지 안띠는지 여부를 검정하는 방법 대해 알아보고, R의 (1) shapiro.test() 함수를 이용한 Shapiro-Wilk 검정, (2) Histogram, Kernel Density Plot 그래프를 이용한 정규성 확인, (3) Q-Q plot 그래프를 이용한 정규성 확인 (qqnorm(), qqline() 함수)를 사용해 예를 들어보겠습니다.

실습할 데이터셋으로는 UsingR 패키지에 들어있는 cfb 데이터 프레임의 소득(income)을 가지고 예를 들어보겠습니다. cfb 데이터셋은 소비자 재정에 관한 설문조사 샘플 데이터로서, 14개의 변수와 1000명의 관측치가 들어있습니다.

> installed.packages("UsingR")
     Package LibPath Version Priority Depends Imports LinkingTo Suggests Enhances License
     License_is_FOSS License_restricts_use OS_type Archs MD5sum NeedsCompilation Built
> library(UsingR)
필요한 패키지를 로딩중입니다: MASS
필요한 패키지를 로딩중입니다: HistData
필요한 패키지를 로딩중입니다: Hmisc
필요한 패키지를 로딩중입니다: grid
필요한 패키지를 로딩중입니다: lattice
필요한 패키지를 로딩중입니다: survival
필요한 패키지를 로딩중입니다: Formula
필요한 패키지를 로딩중입니다: ggplot2
Find out what's changed in ggplot2 with
news(Version == "1.0.1", package = "ggplot2")

다음의 패키지를 부착합니다: ‘Hmisc’

The following objects are masked from ‘package:base’:

    format.pval, round.POSIXt, trunc.POSIXt, units


다음의 패키지를 부착합니다: ‘UsingR’

The following object is masked from ‘package:ggplot2’:

    movies

The following object is masked from ‘package:survival’:

    cancer

> str(cfb)
'data.frame':	1000 obs. of  14 variables:
 $ WGT     : num  5750 5871 8044 6093 7162 ...
 $ AGE     : num  54 40 35 55 40 82 26 50 71 70 ...
 $ EDUC    : num  14 12 14 12 12 12 16 14 12 6 ...
 $ INCOME  : num  66814 42144 25698 35977 39061 ...
 $ CHECKING: num  6000 400 1000 2600 1000 1000 3000 3100 1000 50 ...
 $ SAVING  : num  2000 0 160 19100 8300 0 0 0 0 0 ...
 $ NMMF    : num  0 0 0 0 0 50000 0 0 0 0 ...
 $ STOCKS  : num  500 0 0 0 3500 0 0 0 0 0 ...
 $ FIN     : num  39600 5400 15460 54700 12800 ...
 $ VEHIC   : num  6400 21000 2000 18250 9100 ...
 $ HOMEEQ  : num  84000 8000 12000 90000 47000 175000 0 22000 15000 0 ...
 $ OTHNFIN : num  0 0 0 0 0 0 0 0 0 0 ...
 $ DEBT    : num  40200 58640 19610 8000 21000 ...
 $ NETWORTH: num  170800 17760 9850 284950 268900 ...

(1) Shapiro-Wilk normality test : shapiro.test(x)

귀무가설 H0 : 모집단은 정규분포를 따른다

대립가설 H1 : 모집단은 정규분포를 따르지 않는다

[ Shapiro-Wilk test W statistics ]

- source : https://en.wikipedia.org/wiki/Shapiro%E2%80%93Wilk_test

> # Shapiro-Wilk normality test
> shapiro.test(cfb$INCOME)

	Shapiro-Wilk normality test

data:  cfb$INCOME
W = 0.36883, p-value < 2.2e-16

P-value 값이 매우 작으므로 귀무가설 (H0 : 수입은 정규분포를 따른다)을 기각하고, 대립가설 (H1: 수입은 정규분포를 따르지 않는다)를 채택하게 되겠습니다.

통계량만을 가지고 정규성 여부를 검정하는 것은 부족하며, 반드시 그래프를 통해서 정규성 여부를 병행해서 확인할 필요가 있습니다. Histogram, Kernel density plot 이 빠르고 쉽게, 직관적으로 단변량 분포를 확인할 수 있는 방법입니다.

(2) Histogram and Kernel Density Plot 을 이용한 정규성 확인 : hist(x), lines(density(x))

> # Histogram : hist()
> hist(cfb$INCOME, breaks=100)

아래의 Kernel Density Plot에서는 Y축이 Frequency가 아니라 Probability 임에 유의하시기 바랍니다. (hist(freq=FALSE)) 옵션 사용)

> # Kernel Density Plot : hist(freq=FALSE), lines(density())
> hist(cfb$INCOME, freq=FALSE, breaks=100, main="Kernel Density Plot of cfb$INCOME")
> lines(density(cfb$INCOME), col="blue", lwd=3)

위의 두 개의 그래프를 보면 수입(INCOME)은 F-분포(F-distribution) 혹은 멱함수 분포 (Power-law distribution)을 따르고 있음을 알 수 있습니다. 왼쪽으로 심하게 치우쳐 있고, 오른쪽으로 꼬리가 긴 분포이므로 정규분포와는 거리가 멀다고 하겠습니다.

(3) 분위수-분위수 그림 Q-Q plot을 이용한 정규성 확인 : qqnorm(x), qqline(x) 함수

> # Q-Q plot : qqnorm(), qqline()
> qqnorm(cfb$INCOME)
> qqline(cfb$INCOME)

분위수-분위수 그램 (Q-Q plot : Quantile-Quantile plot) 은 (q(i), X(i))를 2차원 평면에 산점도를 그린 것입니다.

[ Empirical Quantile vs. Theoretical Quantile ]

위의 Q-Q plot 해석 방법은, 경험분포(empirical distribution)와 이론적 분포(theoretical distribution)가 서로 근접하게 분포하고 있으면 정규성을 띤다고 평가하며, 반대로 경험분포와 이론적 분포가 서로 멀어질 수록 정규분포를 띠지 않는다고 평가합니다. 다시말하면, 정규확률 산점도 그림에서 대각선(qqline())을 기준으로 산점도 점들이 가깝게 선형을 이루며 붙어 있으면 정규성을 띤다고 평가하고, 그렇지 않으면 정규성을 띠지 않는다고 보는 것입니다.

R을 사용해서 한쪽으로 치우친 데이터에 대한 정규분포로의 변환 방법에 대해서는 R 데이터 변환 (2) 정규분포화 log(), sqrt() 을 참고하시기 바랍니다.

모집단의 분포 형태에 따른 대략적인 정규분포 변환 방법은 아래 표와 같습니다

[ 분포 형태별 정규분포 변환 방법 ]

distribution before transformation	transformation function	distribution after transformation
left	X^3	normal distribution (bell shape)
mild left	X^2
mild right	sqrt(X)
right	ln(X)
severe right	1/X

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

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

728x90

저작자표시 비영리 변경금지

'R 분석과 프로그래밍 > R 통계분석' 카테고리의 다른 글

R 이변량 범주형 자료의 정리와 요약 - 분할표(contingency table) : table(), xtabs(), prop.table(), margin.table(), addmargins(), CrossTable() (0)	2015.10.10
R 일변량 범주형 자료의 정리와 요약 - 도수분포표(frequency distribution table) : table(), xtabs() (2)	2015.10.10
R 단일 모집단 중심에 대한 비모수 검정 : wilcox.test() (0)	2015.10.04
R 단일 모집단의 모비율에 대한 신뢰구간 추정과 검정(population proportion test) : prop.test() (0)	2015.10.03
R 단일 모집단의 모분산(population variance)에 대한 신뢰구간 추정과 검정 : chisq test (3)	2015.10.03

Posted by Rfriend

이전 1 ··· 17 18 19 20 21 22 23 ··· 29 다음

R, Python 분석과 프로그래밍의 친구 (by R Friend)

'R 분석과 프로그래밍'에 해당되는 글 289건

R (1) 독립된 두 표본의 모평균 차이에 대한 추정과 검정 : t.test()

'R 분석과 프로그래밍 > R 통계분석' 카테고리의 다른 글

R 두 연속형 변수 상관관계 : 공분산(covariance) cov(), 상관계수(correlation coefficient) cor()

'R 분석과 프로그래밍 > R 통계분석' 카테고리의 다른 글

R 연속형 변수 (2) 그룹별(요인별) 요약통계 비교하기 : tapply(var, factor, summary), by(), aggregate(), summaryBy(), describeBy()

'R 분석과 프로그래밍 > R 통계분석' 카테고리의 다른 글

R 연속형 변수 (1) 요약통계 한번에 보기 : summary(), stat.desc(), describe()

'R 분석과 프로그래밍 > R 통계분석' 카테고리의 다른 글

R 일변량 연속형 자료 요약 - (3) 분포형태와 대칭정도(distrubution, symmetry) : 왜도(skewness), 첨도(kurtosis), 분위수-분위수(Quantile-Qunantile)

'R 분석과 프로그래밍 > R 통계분석' 카테고리의 다른 글

R 일변량 연속형 자료 요약 - (2) 퍼짐 정도(dispersion) : 분산(variance), 표준편차(sd), 변이계수, 범위(range), IQR, 백분위수(percentile)

'R 분석과 프로그래밍 > R 통계분석' 카테고리의 다른 글

R 일변량 연속형 자료 요약 - (1) 중심화 경향 : 산술평균(mean), 중앙값(median), 최빈값, 기하평균, CAGR, 조화평균, 가중평균

'R 분석과 프로그래밍 > R 통계분석' 카테고리의 다른 글

R 이변량 범주형 자료의 정리와 요약 - 분할표(contingency table) : table(), xtabs(), prop.table(), margin.table(), addmargins(), CrossTable()

'R 분석과 프로그래밍 > R 통계분석' 카테고리의 다른 글

R 일변량 범주형 자료의 정리와 요약 - 도수분포표(frequency distribution table) : table(), xtabs()

'R 분석과 프로그래밍 > R 통계분석' 카테고리의 다른 글

R 단일 모집단 분포의 정규성 검정 : shapiro.test(), qqnorm(), qqline()

'R 분석과 프로그래밍 > R 통계분석' 카테고리의 다른 글

카테고리

태그목록

티스토리툴바