데이터분석/R

[R 통계분석] 기술통계학

버섯도리 2022. 8. 9. 06:09

cafedata.csv
0.00MB

> ## 00) data 준비

> # 통계청 마이크로데이터 통합서비스(https://mdis.kostat.go.kr/) - 인구 > 인구주택총조사 > 1%_인구사항(제공) > 2010
> data <- read.csv("./data/ch02.csv", na.strings = c("."))
> str(data)
'data.frame': 468284 obs. of  6 variables:
 $ 성별코드         : int  1 1 1 1 1 1 1 1 1 1 ...
 $ 만연령           : int  0 0 0 0 0 0 0 0 0 0 ...
 $ 가구주관계코드   : int  3 3 3 3 3 3 3 3 3 3 ...
 $ 교육정도코드     : int  1 1 1 1 1 1 1 1 1 1 ...
 $ 총출생아수_남아수: int  NA NA NA NA NA NA NA NA NA NA ...
 $ 총출생아수_여아수: int  NA NA NA NA NA NA NA NA NA NA ...
> unique(data$성별코드)
[1] 1 2
> data$성별코드 <- factor(data$성별코드, levels = c(1, 2),
+                     labels = c("남자","여자"))
> sort(unique(data$가구주관계코드))
 [1]  1  2  3  4  5  6  7  8  9 10 11 12 13 14
> data$가구주관계코드 <- factor(data$가구주관계코드, levels=1:14,
+                        labels = c("가구주", "가구주의 배우자", "자녀",
+                                   "자녀의 배우자", "가구주의 부모",
+                                   "배우자의 부모", "손자녀, 그 배우자",
+                                   "증손자녀, 그 배우자", "조부모",
+                                   "형제자매, 그 배우자",
+                                   "형제자매의 자녀, 그 배우자",
+                                   "부모의 형제자매, 그 배우자", "기타 친인척",
+                                   "그외같이사는사람"))
> sort(unique(data$교육정도코드))
[1] 1 2 3 4 5 6 7 8
> data$교육정도코드 <- factor(data$교육정도코드, levels=1:8,
+                       labels = c("안 받았음", "초등학교", "중학교",
+                                  "고등학교", "대학-4년제 미만", "대학-4년제 이상",
+                                  "석사과정", "박사과정"))
> str(data)
'data.frame': 468284 obs. of  6 variables:
 $ 성별코드         : Factor w/ 2 levels "남자","여자": 1 1 1 1 1 1 1 1 1 1 ...
 $ 만연령           : int  0 0 0 0 0 0 0 0 0 0 ...
 $ 가구주관계코드   : Factor w/ 14 levels "가구주","가구주의 배우자",..: 3 3 3 3 3 3 3 3 3 3 ...
 $ 교육정도코드     : Factor w/ 8 levels "안 받았음","초등학교",..: 1 1 1 1 1 1 1 1 1 1 ...
 $ 총출생아수_남아수: int  NA NA NA NA NA NA NA NA NA NA ...
 $ 총출생아수_여아수: int  NA NA NA NA NA NA NA NA NA NA ...
> colnames(data)[c(1,3,4)]
[1] "성별코드"       "가구주관계코드" "교육정도코드"  
> colnames(data)[c(1,3,4)] <- c("성별", "가구주관계", "교육정도")
> str(data)
'data.frame': 468284 obs. of  6 variables:
 $ 성별             : Factor w/ 2 levels "남자","여자": 1 1 1 1 1 1 1 1 1 1 ...
 $ 만연령           : int  0 0 0 0 0 0 0 0 0 0 ...
 $ 가구주관계       : Factor w/ 14 levels "가구주","가구주의 배우자",..: 3 3 3 3 3 3 3 3 3 3 ...
 $ 교육정도         : Factor w/ 8 levels "안 받았음","초등학교",..: 1 1 1 1 1 1 1 1 1 1 ...
 $ 총출생아수_남아수: int  NA NA NA NA NA NA NA NA NA NA ...
 $ 총출생아수_여아수: int  NA NA NA NA NA NA NA NA NA NA ...
> save.image("data.rda")
> summary(data)
   성별            만연령                    가구주관계                교육정도      총출생아수_남아수 총출생아수_여아수
 남자:226965   Min.   : 0.00   가구주             :179293   고등학교       :134246   Min.   : 0.00     Min.   : 0.00    
 여자:241319   1st Qu.:22.00   자녀               :145533   대학-4년제 이상: 81110   1st Qu.: 1.00     1st Qu.: 0.00    
               Median :40.00   가구주의 배우자    :106311   초등학교       : 80710   Median : 1.00     Median : 1.00    
               Mean   :39.34   가구주의 부모      : 12069   중학교         : 55704   Mean   : 1.32     Mean   : 1.23    
               3rd Qu.:55.00   손자녀, 그 배우자  :  7832   안 받았음      : 51085   3rd Qu.: 2.00     3rd Qu.: 2.00    
               Max.   :85.00   형제자매, 그 배우자:  5332   대학-4년제 미만: 50753   Max.   :12.00     Max.   :11.00    
                               (Other)            : 11914   (Other)        : 14676   NA's   :310308    NA's   :310308   
> head(subset(data, subset = is.na(총출생아수_남아수) & !is.na(총출생아수_여아수)))
[1] 성별          만연령            가구주관계    교육정도      총출생아수_남아수 총출생아수_여아수
<0 행> <또는 row.names의 길이가 0입니다>
> head(subset(data, subset = !is.na(총출생아수_남아수) & is.na(총출생아수_여아수)))
[1] 성별          만연령            가구주관계    교육정도      총출생아수_남아수 총출생아수_여아수
<0 행> <또는 row.names의 길이가 0입니다>

 

Section 01. 그래프

> ## 01) 두 변수 간의 관계를 나타내는 산점도

> str(cars)
'data.frame': 50 obs. of  2 variables:
 $ speed: num  4 4 7 7 8 9 10 10 10 11 ...
 $ dist : num  2 10 4 22 16 10 18 26 34 17 ...
plot(cars$speed, cars$dist,
+      main = "속도와 제동거리", xlab = "속도(mph)", ylab = "제동거리(ft)",
+      pch = 1, col = "red")


> Nile
Time Series:
Start = 1871 
End = 1970 
Frequency = 1 
  [1] 1120 1160  963 1210 1160 1160  813 1230 1370 1140  995  935 1110  994 1020  960 1180  799  958 1140 1100 1210 1150
 [24] 1250 1260 1220 1030 1100  774  840  874  694  940  833  701  916  692 1020 1050  969  831  726  456  824  702 1120
 [47] 1100  832  764  821  768  845  864  862  698  845  744  796 1040  759  781  865  845  944  984  897  822 1010  771
 [70]  676  649  846  812  742  801 1040  860  874  848  890  744  749  838 1050  918  986  797  923  975  815 1020  906
 [93]  901 1170  912  746  919  718  714  740
> str(Nile)
 Time-Series [1:100] from 1871 to 1970: 1120 1160 963 1210 1160 1160 813 1230 1370 1140 ...
plot(Nile,
+      main = "Nile강의 연도별 유량 변화", xlab = "연도", ylab = "유량")

plot(Nile, type = "p",
+      main = "Nile강의 연도별 유량 변화", xlab = "연도", ylab = "유량")

>

> ## 02) 막대그래프와 히스토그램

> load("data.rda")
> tableBabyM <- table(data$총출생아수_남아수)
> tableBabyM

    0     1     2     3     4     5     6     7     8     9    10    11    12 
30788 69624 41010 11165  3667  1228   346   104    21     8     4    10     1 
barplot(tableBabyM,
+         main = "출생아(남자)별 빈도", xlab = "출생아수", ylab = "빈도")


> tableSexEdu <- table(data$성별, data$교육정도)
> tableSexEdu
      
       안 받았음 초등학교 중학교 고등학교 대학-4년제 미만 대학-4년제 이상 석사과정 박사과정
  남자     19161    34214  26588    66548           25673           45530     7107     2144
  여자     31924    46496  29116    67698           25080           35580     4634      791
barplot(tableSexEdu, legend.text = T, col = c("orange", "green"),
+         main = "학력에 따른 성별 인원수", xlab = "학력", ylab = "빈도")

barplot(tableSexEdu, legend.text = T, col = c("orange", "green"), beside = T,
+         main = "학력에 따른 성별 인원수", xlab = "학력", ylab = "빈도")


hist(data$만연령, main = "연령별 분포", xlab = "연령", ylab = "빈도")

hist(data$만연령, breaks = c(seq(0, 90, 10)), right = F,
+      main = "연령별 분포", xlab = "연령", ylab = "빈도")

>

> ## 03) 원도표

> load("data.rda")
> tableEdt <- table(data$교육정도)
> tableEdt  

      안 받았음        초등학교          중학교        고등학교 대학-4년제 미만 대학-4년제 이상        석사과정 
          51085           80710           55704          134246           50753           81110           11741 
       박사과정 
           2935 
pie(tableEdt, main = "학력수준별 비중")

>

 

Section 02. 모수와 통계량


> ranicafe <- read.csv("./data/cafedata.csv", header=T, na.strings="na", stringsAsFactors=FALSE )
> ranicafe <- na.omit(ranicafe)

> str(ranicafe)
'data.frame': 47 obs. of  22 variables:
 $ t                        : int  1 2 3 4 5 6 7 8 9 10 ...
 $ Date                     : chr  "2010-01-19" "2010-01-20" "2010-01-21" "2010-01-22" ...
 $ Day.Code                 : int  2 3 4 5 1 2 3 4 5 1 ...
 $ Day.of.Week              : chr  "Tue" "Wed" "Thu" "Fri" ...
 $ Bread.Sand.Sold          : int  5 6 8 4 3 7 6 0 3 2 ...
 $ Bread.Sand.Waste         : int  3 8 2 2 0 1 6 0 4 6 ...
 $ Wraps.Sold               : int  25 7 14 5 10 5 19 7 4 13 ...
 $ Wraps.Waste              : int  5 17 0 7 0 3 3 0 9 3 ...
 $ Muffins.Sold             : int  5 3 4 5 8 1 6 6 0 3 ...
 $ Muffins.Waste            : int  1 5 0 0 0 0 0 1 4 0 ...
 $ Cookies.Sold             : int  5 1 1 3 3 5 10 0 3 6 ...
 $ Cookies.Waste            : int  3 6 0 1 0 0 0 0 2 0 ...
 $ Fruit.Cup.Sold           : int  1 0 0 3 2 2 2 0 1 2 ...
 $ Fruit.Cup.Waste          : int  4 3 3 0 0 0 0 0 1 0 ...
 $ Chips                    : int  12 0 0 20 0 4 2 20 3 16 ...
 $ Juices                   : int  8 0 13 0 5 4 5 6 4 7 ...
 $ Sodas                    : int  20 13 23 13 13 33 15 27 12 19 ...
 $ Coffees                  : int  41 33 34 27 20 23 32 31 30 27 ...
 $ Total.Soda.and.Coffee    : int  61 46 57 40 33 56 47 58 42 46 ...
 $ Sales                    : num  200 196 103 163 102 ...
 $ Max.Daily.Temperature..F.: int  36 34 39 40 36 26 34 33 20 37 ...
 $ Total.Items.Wasted       : int  16 39 5 10 0 4 9 1 20 9 ...
 - attr(*, "na.action")= 'omit' Named int 34
  ..- attr(*, "names")= chr "34"
> head(ranicafe)
  t       Date Day.Code Day.of.Week Bread.Sand.Sold Bread.Sand.Waste Wraps.Sold Wraps.Waste Muffins.Sold Muffins.Waste
1 1 2010-01-19        2         Tue               5                3         25           5            5             1
2 2 2010-01-20        3         Wed               6                8          7          17            3             5
3 3 2010-01-21        4         Thu               8                2         14           0            4             0
4 4 2010-01-22        5         Fri               4                2          5           7            5             0
5 5 2010-01-25        1         Mon               3                0         10           0            8             0
6 6 2010-01-26        2         Tue               7                1          5           3            1             0
  Cookies.Sold Cookies.Waste Fruit.Cup.Sold Fruit.Cup.Waste Chips Juices Sodas Coffees Total.Soda.and.Coffee  Sales
1            5             3              1               4    12      8    20      41                    61 199.95
2            1             6              0               3     0      0    13      33                    46 195.74
3            1             0              0               3     0     13    23      34                    57 102.68
4            3             1              3               0    20      0    13      27                    40 162.88
5            3             0              2               0     0      5    13      20                    33 101.76
6            5             0              2               0     4      4    33      23                    56 186.94
  Max.Daily.Temperature..F. Total.Items.Wasted
1                        36                 16
2                        34                 39
3                        39                  5
4                        40                 10
5                        36                  0
6                        26                  4

> library(ggplot2)
> ggplot(ranicafe, aes(Coffees)) + 
+   geom_bar(fill="gray") + 
+   ggtitle("라니의 카페 커피 판매량") + 
+   theme(plot.title = element_text(size = 20, face="bold")) +
+   xlim(0, 50) + xlab("판매량") +
+   ylab("횟수") + scale_y_continuous(breaks=0:10)


> ## 04) 최대값과 최소값

> sort(ranicafe$Coffees)[1]
[1] 3
> sort(ranicafe$Coffees, decreasing = T)[1]
[1] 48
> min(ranicafe$Coffees)
[1] 3
> max(ranicafe$Coffees)
[1] 48

> # 최빈값
> rc <- ranicafe$Coffees
> stem(rc)

  The decimal point is 1 digit(s) to the right of the |

  0 | 34444
  0 | 5688
  1 | 01134
  1 | 668
  2 | 001123344
  2 | 55677789
  3 | 001112334
  3 | 55
  4 | 1
  4 | 8

> #-> 4가 4번 등장하여 최빈값은 4
> tableRc <- table(rc)
> names(tableRc)[which(tableRc == max(tableRc))]
[1] "4"

> ## 05) 평균
> weight <- 1 / length(rc)
> sum(rc * weight)
[1] 21.51064
> mean(rc)
[1] 21.51064

> ## 07) 중앙값
> (median.idx <- (length(rc)+1) / 2)
[1] 24
> (rc.srt <- sort(rc))
 [1]  3  4  4  4  4  5  6  8  8 10 11 11 13 14 16 16 18 20 20 21 21 22 23 23 24 24 25 25 26 27 27 27 28 29 30 30 31 31 31
[40] 32 33 33 34 35 35 41 48
> rc.srt[median.idx]
[1] 23
> median(rc)
[1] 23

> ## 09) 개별 관찰값과 평균과의 차이에 대한 평균
> height <- c(164, 166, 168, 170, 172, 174, 176)
> (height.m <- mean(height))
[1] 170
> (height.dev <- height - height.m)
[1] -6 -4 -2  0  2  4  6
> sum(height.dev)
[1] 0

> ## 10) 편차 제곱의 평균 구하기
> (height.dev2 <- height.dev ^ 2)
[1] 36 16  4  0  4 16 36
> sum(height.dev2)
[1] 112
> mean(height.dev2)
[1] 16

> ## 11) 표준편차 구하기
> sqrt(mean(height.dev2))
[1] 4

> ## 12) 분산과 표준편차 구하기
> var(height)
[1] 18.66667
> sd(height)
[1] 4.320494
> # R은 대부분 표본조사라고 판단하고 갯수(n)가 아닌 자유도(n-1)로 나눈다.

> ## 13) '라나의 카페' 커피 판매량에 대한 평균과 표준편차
> rc.m <- mean(rc)
> rc.sd <- sd(rc)
> cat("커피 판매량", round(rc.m, 1), "±", round(rc.sd, 2), "잔")
커피 판매량 21.5 ± 11.08 잔

> # 정규분포
> par(mar=c(2, 0, 0, 0))
> x <- seq(-3, 3, by=0.001)
> plot(x, dnorm(x), type="l", lwd=2, xlab="", ylab="", axes=F)
> axis(1)
> lines(c(0, 0), c(0, dnorm(0)), lty=1, lwd=1.5)
> lines(c(-1, -1), c(0, dnorm(-1)), lty=3, lwd=2)
> lines(c(1, 1), c(0, dnorm(1)), lty=3, lwd=2)


> # 표준편차가 다른 두 자료
> par(mar=c(2, 0, 0, 0))
> x <- seq(158, 182, by=0.01)
> plot(x, dnorm(x, mean=170, sd=4), type="l", lwd=2, ylim=c(0, 0.2), axes=F)
> lines(x, dnorm(x, mean=170, sd=2), lty=2, lwd=2)
> axis(1)
> legend("topright",
+        lwd=2,
+        lty=c(1,2),
+        c(paste("sd :", 4), paste("sd :", 2))
+ )

>

> ## 14) 변동계수 구하기
> rj <- ranicafe$Juices
> (rj.m <- mean(rj))
[1] 4.93617
> (rj.sd <- sd(rj))
[1] 3.703138
> (rc.cv <- round(rc.sd / rc.m, 3))
[1] 0.515
> (rj.cv <- round(rj.sd / rj.m, 3))
[1] 0.75
>

> ## 15) 사분위수 범위와 상자도표
> (qs <- quantile(rc))
  0%  25%  50%  75% 100% 
   3   12   23   30   48 
> (qs[4] - qs[2])
75% 
 18 
IQR(rc)
[1] 18
> bp <- boxplot(rc, main = "커피 판매량에 대한 상자도표")
> axis(2)
> md <- median(rc)
> arrows(1.4, qs[4], 1.22, qs[4], angle=30, length=0.1)
> text(1.42, qs[4], "Q3", adj=0)
> arrows(1.4, md, 1.22, md, angle=30, length=0.1)
> text(1.42, md, "Q2", adj=0)
> arrows(1.4, qs[2], 1.22, qs[2], angle=30, length=0.1)
> text(1.42, qs[2], "Q1", adj=0)


> # 제동거리의 히스토그램
> hist(cars$dist, breaks = seq(0, 120, 10))


> ## 16) 이상치 판별
> ( Q <- quantile(cars$dist))
  0%  25%  50%  75% 100% 
   2   26   36   56  120 
> (ll <- Q[2] - 1.5*IQR(cars$dist))
25% 
-19 
> (ul <- Q[4] + 1.5*IQR(cars$dist))
75% 
101 
> cars$dist[cars$dist < ll]
numeric(0)
> cars$dist[cars$dist > ul]
[1] 120
boxplot(cars$dist, main = "Boxplot of Distance")

 

 

 

 

 

 

 

출처 : 이윤환, ⌜제대로 알고 쓰는 R 통계분석⌟, 한빛아카데미, 2020