1. 추론통계
1.1 정규분포
x=seq(-5,5,length=1000)
y=dnorm(x,mean = 0,sd=1)
plot(x,y,type="l",lwd=1)
den.norm=function(x)dnorm(x,mean=0,sd=2)
den.norm2=function(x)dnorm(x,mean=0,sd=3)
curve(den.norm,col="red",add=T,lty=2)
curve(den.norm2,col="blue",add=T,lty=2)1.2 표본 추출
1.2.1 단순 임의 추출 1 - 확률이 같은 경우
x=c(1:10)
sample(x,5) # 비복원 추출 ## [1] 8 9 10 7 4
sample(x,5,replace = TRUE) # 복원 추출## [1] 8 2 3 5 2
1.2.2 단순 임의 추출 2 - 확률이 다른 경우
x=c("A","B","C","D","E")
prob=c(6,1,5,1,3) # 가중치 부여
sample(x,3,replace = TRUE,prob = prob) ## [1] "A" "C" "C"
1.2.3 층화 임의 추출 - strata()
library(sampling)sampledata=strata(c("Species"),size=c(2,3,4),method="srswor",data=iris)
getdata(iris,sampledata)## Sepal.Length Sepal.Width Petal.Length Petal.Width Species ID_unit
## 14 4.3 3.0 1.1 0.1 setosa 14
## 29 5.2 3.4 1.4 0.2 setosa 29
## 66 6.7 3.1 4.4 1.4 versicolor 66
## 73 6.3 2.5 4.9 1.5 versicolor 73
## 91 5.5 2.6 4.4 1.2 versicolor 91
## 103 7.1 3.0 5.9 2.1 virginica 103
## 120 6.0 2.2 5.0 1.5 virginica 120
## 128 6.1 3.0 4.9 1.8 virginica 128
## 132 7.9 3.8 6.4 2.0 virginica 132
## Prob Stratum
## 14 0.04 1
## 29 0.04 1
## 66 0.06 2
## 73 0.06 2
## 91 0.06 2
## 103 0.08 3
## 120 0.08 3
## 128 0.08 3
## 132 0.08 3
class(sampledata)## [1] "data.frame"
1.2.4 연습문제 1
item=c("A","B","C","D","E")
sales=c(30,5,25,5,15)
df=data.frame(item,sales)
df## item sales
## 1 A 30
## 2 B 5
## 3 C 25
## 4 D 5
## 5 E 15
df$sales/sum(df$sales)## [1] 0.3750 0.0625 0.3125 0.0625 0.1875
dfsample=sample(df$item,replace = TRUE,size = 3,prob = df$sales/sum(df$sales))
dfsample## [1] D C A
## Levels: A B C D E
1.2.5 타이타닉 연습문제 1
- 샘플 사이즈 5%,10%,50% 일때 생존률 비교 프로그램 작성 (for문 이용)
titanic=read.csv("data/Titanic/train.csv")data=titanic$Survived
size=length(data)
weight=c(0.05,0.1,0.5)
wsize=c(size*weight)
for(i in 1:length(weight))
{
survive=sample(data,size=wsize[i],replace=FALSE)
print(paste("추출 비율이",weight[i]*100,"%","일때,"))
print(paste("표본 생존률은",(round(sum(survive)/length(survive)*100)),"%","입니다."))
}## [1] "<U+CD94><U+CD9C> <U+BE44><U+C728><U+C774> 5 % <U+C77C><U+B54C>,"
## [1] "<U+D45C><U+BCF8> <U+C0DD><U+C874><U+B960><U+C740> 39 % <U+C785><U+B2C8><U+B2E4>."
## [1] "<U+CD94><U+CD9C> <U+BE44><U+C728><U+C774> 10 % <U+C77C><U+B54C>,"
## [1] "<U+D45C><U+BCF8> <U+C0DD><U+C874><U+B960><U+C740> 43 % <U+C785><U+B2C8><U+B2E4>."
## [1] "<U+CD94><U+CD9C> <U+BE44><U+C728><U+C774> 50 % <U+C77C><U+B54C>,"
## [1] "<U+D45C><U+BCF8> <U+C0DD><U+C874><U+B960><U+C740> 39 % <U+C785><U+B2C8><U+B2E4>."
1.2.6 타이타닉 연습문제 2
- 남성과 여성 생존률 평균
- 남성과 여성 층화 임의 추출 50명씩 추출, 평균 생존률 (비복원)
- 남성, 여성 데이터 생성
titanicmale=titanic[titanic$Sex=="male",] # 타이타닉 남성 데이터 생성
titanicfemale=titanic[titanic$Sex=="female",] # 타이타닉 여성 데이터 생성- 모집단 남성 생존률
malesurvive=sum(titanicmale$Survive)/length(titanicmale$Survive) # 타이타닉 남성 생존률
malesurvive## [1] 0.1889081
- 모집단 여성 생존률
femalesurvive=sum(titanicfemale$Survive)/length(titanicfemale$Survive) # 타이타닉 여성 생존률
femalesurvive## [1] 0.7420382
- 층화 임의 추출
stitanic=strata(c("Sex"),size = c(50,50),method = "srswor",data=titanic) # 층화 임의 추출
head(stitanic)## Sex ID_unit Prob Stratum
## 1 male 1 0.08665511 1
## 6 male 6 0.08665511 1
## 46 male 46 0.08665511 1
## 49 male 49 0.08665511 1
## 52 male 52 0.08665511 1
## 70 male 70 0.08665511 1
stitanic=getdata(titanic,stitanic)- 층화 임의 추출 후 남성/여성 데이터 생성
stitanicmale=stitanic[stitanic$Sex=="male",]
stitanicfemale=stitanic[stitanic$Sex=="female",]- 표본 남성 생존률
smalesurvive=sum(stitanicmale$Survive)/length(stitanicmale$Survive) # 표본 남성 생존률
smalesurvive## [1] 0.2
- 표본 여성 생존률
sfemalesurvive=sum(stitanicfemale$Survive)/length(stitanicfemale$Survive) # 표본 여성 생존률
sfemalesurvive## [1] 0.74
1.3 가설검정
- 데이터 불러오기
library(MASS)
data("survey")
head(survey)## Sex Wr.Hnd NW.Hnd W.Hnd Fold Pulse Clap Exer Smoke Height
## 1 Female 18.5 18.0 Right R on L 92 Left Some Never 173.00
## 2 Male 19.5 20.5 Left R on L 104 Left None Regul 177.80
## 3 Male 18.0 13.3 Right L on R 87 Neither None Occas NA
## 4 Male 18.8 18.9 Right R on L NA Neither None Never 160.00
## 5 Male 20.0 20.0 Right Neither 35 Right Some Never 165.00
## 6 Female 18.0 17.7 Right L on R 64 Right Some Never 172.72
## M.I Age
## 1 Metric 18.250
## 2 Imperial 17.583
## 3 <NA> 16.917
## 4 Metric 20.333
## 5 Metric 23.667
## 6 Imperial 21.000
1.3.1 독립성 검정 - chisq.test()
- 데이터 테이블 도수분포표 작성
xt=xtabs(~Sex+Exer,data=survey) # 데이터테이블의 도수분포표 작성
xt## Exer
## Sex Freq None Some
## Female 49 11 58
## Male 65 13 40
- 카이제곱 테스트
chisq.test(xt) # 카이제곱 테스트로 독립성 검정 ##
## Pearson's Chi-squared test
##
## data: xt
## X-squared = 5.7184, df = 2, p-value = 0.05731
1.3.2 피셔의 정확 검정 - fisher.test()
- 데이터가 한쪽에 몰려있을 때 사용
- 카이제곱 테스트가 정확하지 않을 수 있다는 메시지 출력
xt=xtabs(~W.Hnd+Clap,data=survey)
chisq.test(xt) ## Warning in chisq.test(xt): Chi-squared approximation may be incorrect
##
## Pearson's Chi-squared test
##
## data: xt
## X-squared = 19.252, df = 2, p-value = 6.598e-05
fisher.test(xt)##
## Fisher's Exact Test for Count Data
##
## data: xt
## p-value = 0.0001413
## alternative hypothesis: two.sided
1.3.3 mcnemar 검정 - mcnemar.test()
- 어떤 사건의 전후로 응답자의 성향등을 조사 할 때 사용
- 예시로는 TV 유세 후 지지율 변화
- 데이터는 매트릭스 타입이어야 함
performance=matrix(c(794,86,150,570),nrow=2,
dimnames=list("1st survey"=c("approve","disapprove"),
"2nd survet"=c("approve","disapprove")))
performance## 2nd survet
## 1st survey approve disapprove
## approve 794 150
## disapprove 86 570
mcnemar.test(performance)##
## McNemar's Chi-squared test with continuity correction
##
## data: performance
## McNemar's chi-squared = 16.818, df = 1, p-value = 4.115e-05
1.3.4 카이제곱 테스트 - chisq.test(table,p)
- 데이터가 특정 분포를 따르는지 살펴볼 때 사용
- 결과 해석 : p-value가 작으므로 주어진 분포(왼손과 오른손의 비율이 0.3/0.7)를 따르지 않는다.
table=table(survey$W.Hnd)
chisq.test(table,p=c(0.3,0.7))##
## Chi-squared test for given probabilities
##
## data: table
## X-squared = 56.252, df = 1, p-value = 6.376e-14
1.3.5 샤피로 윌크 검정 - shapiro.test()
- 표본이 정규 분포로부터 추출 된 것인지 테스트
- p-값이 크므로 귀무가설 채택 -> 주어진 표본은 정규분포를 따른다.
x=rnorm(1000)
shapiro.test(x)##
## Shapiro-Wilk normality test
##
## data: x
## W = 0.99911, p-value = 0.9238
1.3.6 콜모고로프 스미르노프 검정 - ks.test()
- 두 데이터가 차이가 있는지 검정
- 귀무가설 : 같은 분포
- p-value가 크므로 귀무가설 채택 -> 같은 분포를 따른다.
ks.test(rnorm(100),rnorm(90))##
## Two-sample Kolmogorov-Smirnov test
##
## data: rnorm(100) and rnorm(90)
## D = 0.11444, p-value = 0.5202
## alternative hypothesis: two-sided
1.3.7 독립성 검정 연습문제
- 생존여부를 factor로 변환
- chisq 이용하여 성별과 생존여부의 p-value 구하기
- chisq 이용하여 좌석등급과 생존여부의 p-value 구하기
- chisq 이용하여 이름과 생존여부의 p-value 구하기
- 3개의 케이스 귀무가설 기각 여부 판단
- factor 변환
titanic$Survived=as.factor(titanic$Survived) # factor 변환
str(titanic$Survived)## Factor w/ 2 levels "0","1": 1 2 2 2 1 1 1 1 2 2 ...
- 성별 - 생존여부 / 귀무가설 기각
xt1=xtabs(~Sex+Survived,data=titanic) # 성별 - 생존여부 / 귀무가설 기각
chisq.test(xt1) ##
## Pearson's Chi-squared test with Yates' continuity correction
##
## data: xt1
## X-squared = 260.72, df = 1, p-value < 2.2e-16
- 좌석등급 - 생존여부 / 귀무가설 기각
xt2=xtabs(~Pclass+Survived,data=titanic) # 좌석등급 - 생존여부 / 귀무가설 기각
chisq.test(xt2) ##
## Pearson's Chi-squared test
##
## data: xt2
## X-squared = 102.89, df = 2, p-value < 2.2e-16
- 이름 - 생존여부 / 귀무가설 채택
xt3=xtabs(~Name+Survived,data=titanic) # 이름 - 생존여부 / 귀무가설 채택
chisq.test(xt3) ## Warning in chisq.test(xt3): Chi-squared approximation may be incorrect
##
## Pearson's Chi-squared test
##
## data: xt3
## X-squared = 891, df = 890, p-value = 0.4842
1.3.8 독립성 검정 연습문제 2
- 60명 대상자 A,B 그룹별로 대장암에 걸렸는지 검사한 결과
- 그룹과 대장암 발생 여부 사이의 관련성 검증
- 자료 생성
group=c("A","A","B","B")
cancer=c("1.Yes","2.NO","1.Yes","2.NO")
count=c(2,28,5,25)
dat=data.frame(group,cancer,count)
tab=xtabs(count~group+cancer,data=dat)
tab## cancer
## group 1.Yes 2.NO
## A 2 28
## B 5 25
- 피셔 정확검정
fisher.test(tab)##
## Fisher's Exact Test for Count Data
##
## data: tab
## p-value = 0.4238
## alternative hypothesis: true odds ratio is not equal to 1
## 95 percent confidence interval:
## 0.03189243 2.46223985
## sample estimates:
## odds ratio
## 0.3630951
1.3.8 독립성 검정 연습문제 3
방법 1
example2=matrix(c(59,16,6,80),nrow=2,
dimnames=list("after"=c("present","absent"),
"before"=c("present","absent")))
mcnemar.test(example2)##
## McNemar's Chi-squared test with continuity correction
##
## data: example2
## McNemar's chi-squared = 3.6818, df = 1, p-value = 0.05501
방법2
before.group=c("present","present","absent","absent")
after.group=c("present","absent","present","absent")
count=c(59,6,16,80)
dat=data.frame(before.group,after.group,count)
dat## before.group after.group count
## 1 present present 59
## 2 present absent 6
## 3 absent present 16
## 4 absent absent 80
xt=xtabs(count~before.group+after.group,data=dat)
mcnemar.test(xt)##
## McNemar's Chi-squared test with continuity correction
##
## data: xt
## McNemar's chi-squared = 3.6818, df = 1, p-value = 0.05501
1.3.8 독립성 검정 연습문제 4
shapiro.test(titanic$Pclass)##
## Shapiro-Wilk normality test
##
## data: titanic$Pclass
## W = 0.71833, p-value < 2.2e-16
shapiro.test(titanic$Age)##
## Shapiro-Wilk normality test
##
## data: titanic$Age
## W = 0.98146, p-value = 7.337e-08
shapiro.test(titanic$Fare)##
## Shapiro-Wilk normality test
##
## data: titanic$Fare
## W = 0.52189, p-value < 2.2e-16
1.4 상관 분석
- 두 변수 사이의 관련성을 파악하는 방법
- 대표적으로 피어슨 상관계수가 있음
- 피어슨 상관계수는 -1에서 1사이의 값을 갖는다.
1.4.1 상관계수 - 피어슨
cor(iris$Sepal.Length,iris$Sepal.Width) # 피어슨## [1] -0.1175698
1.4.2 상관계수 - 스피어만
- 순위를 사용해서 구하는 방식
- 비선형 상관관계도 측정가능
kor=c(85,90,87,92,95)
eng=c(88,89,68,84,91)
rank(kor) ; rank(eng)## [1] 1 3 2 4 5
## [1] 3 4 1 2 5
m=matrix(c(kor,eng),ncol=2)
cor(m,method = "spearman")## [,1] [,2]
## [1,] 1.0 0.5
## [2,] 0.5 1.0
1.4.3 상관계수 시각화 1 - symnum()
x=iris[,1:4]
head(x)## Sepal.Length Sepal.Width Petal.Length Petal.Width
## 1 5.1 3.5 1.4 0.2
## 2 4.9 3.0 1.4 0.2
## 3 4.7 3.2 1.3 0.2
## 4 4.6 3.1 1.5 0.2
## 5 5.0 3.6 1.4 0.2
## 6 5.4 3.9 1.7 0.4
cor(x)## Sepal.Length Sepal.Width Petal.Length Petal.Width
## Sepal.Length 1.0000000 -0.1175698 0.8717538 0.8179411
## Sepal.Width -0.1175698 1.0000000 -0.4284401 -0.3661259
## Petal.Length 0.8717538 -0.4284401 1.0000000 0.9628654
## Petal.Width 0.8179411 -0.3661259 0.9628654 1.0000000
symnum(cor(x))## S.L S.W P.L P.W
## Sepal.Length 1
## Sepal.Width 1
## Petal.Length + . 1
## Petal.Width + . B 1
## attr(,"legend")
## [1] 0 ' ' 0.3 '.' 0.6 ',' 0.8 '+' 0.9 '*' 0.95 'B' 1
1.4.4 상관계수 시각화 2 - corrgram()
library(corrgram)
corrgram(iris,upper.panel=panel.conf,cex.labels = 2)1.4.5 상관계수 시각화 3 - corrplot()
library(corrplot)## corrplot 0.84 loaded
corrplot 사용예시
how to use corrplot
method
- method=“shade”
corr=cor(iris[,1:4])
corrplot(corr,method = "shade", addCoef.col = "black")- method=“circle”
corrplot(corr,method = "circle", addCoef.col = "black")- method=“ellipse”
corrplot(corr,method = "ellipse", addCoef.col = "black")- method=“number”
corrplot(corr,method = "number", addCoef.col = "black")- method=“square”
corrplot(corr,method = "square", addCoef.col = "black")- method=“pie”
corrplot(corr,method = "pie", addCoef.col = "black")type
- type=“full”
corrplot(corr,method = "ellipse", addCoef.col = "black",type="full")- type=“upper”
corrplot(corr,method = "ellipse", addCoef.col = "black",type="upper")- type=“lower”
corrplot(corr,method = "ellipse", addCoef.col = "black",type="lower")col
corrplot(corr,method = "ellipse", addCoef.col = "black",type="lower",col="yellow")background
corrplot(corr,method = "ellipse", addCoef.col = "black",type="upper",bg="skyblue")diag
- 자기 자신 빼고 싶으면 diag=FALSE
corrplot(corr,method = "ellipse", addCoef.col = "black",type="full",bg="white",diag=FALSE)1.5 집단 차이 분석
1.5.1 일표본 평균 - t.test()
- 한 집단의 평균 검정
- 귀무가설 : 모집단의 평균은 0이다.
- p-값이 크므로 귀무가설 채택 -> 모집단 평균은 0이다.
양측 검정
x=rnorm(30)
t.test(x,alternative = "two.sided",mu=0)##
## One Sample t-test
##
## data: x
## t = 0.33831, df = 29, p-value = 0.7376
## alternative hypothesis: true mean is not equal to 0
## 95 percent confidence interval:
## -0.2830667 0.3952723
## sample estimates:
## mean of x
## 0.05610281
우단측 검정
t.test(x,alternative = "greater",mu=0) # 우단측##
## One Sample t-test
##
## data: x
## t = 0.33831, df = 29, p-value = 0.3688
## alternative hypothesis: true mean is greater than 0
## 95 percent confidence interval:
## -0.225671 Inf
## sample estimates:
## mean of x
## 0.05610281
좌단측 검정
t.test(x,alternative = "less",mu=0)# 좌단측##
## One Sample t-test
##
## data: x
## t = 0.33831, df = 29, p-value = 0.6312
## alternative hypothesis: true mean is less than 0
## 95 percent confidence interval:
## -Inf 0.3378766
## sample estimates:
## mean of x
## 0.05610281
1.5.2 독립 이표본 평균 - t.test()
- 두 모집단의 평균이 같은지 검정
- 데이터 확인
head(sleep)## extra group ID
## 1 0.7 1 1
## 2 -1.6 1 2
## 3 -0.2 1 3
## 4 -1.2 1 4
## 5 -0.1 1 5
## 6 3.4 1 6
tail(sleep)## extra group ID
## 15 -0.1 2 5
## 16 4.4 2 6
## 17 5.5 2 7
## 18 1.6 2 8
## 19 4.6 2 9
## 20 3.4 2 10
sleep2=sleep[,c(1:2)]
head(sleep2)## extra group
## 1 0.7 1
## 2 -1.6 1
## 3 -0.2 1
## 4 -1.2 1
## 5 -0.1 1
## 6 3.4 1
- 모분산 같은지 확인
- 귀무가설 : 분산의 비가 1이다
- p-value가 크므로 귀무가설 채택 -> 등분산
var.test(extra~group,sleep2)##
## F test to compare two variances
##
## data: extra by group
## F = 0.79834, num df = 9, denom df = 9, p-value = 0.7427
## alternative hypothesis: true ratio of variances is not equal to 1
## 95 percent confidence interval:
## 0.198297 3.214123
## sample estimates:
## ratio of variances
## 0.7983426
- t.test (등분산이면 var.equal=TRUE)
- 귀무가설 : 모집단의 평균은 0이다.
- p-값이 크므로 귀무가설 채택 -> 모집단 평균은 0이다.
t.test(extra~group,data=sleep2,paired=FALSE,var.equal=TRUE)##
## Two Sample t-test
##
## data: extra by group
## t = -1.8608, df = 18, p-value = 0.07919
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -3.363874 0.203874
## sample estimates:
## mean in group 1 mean in group 2
## 0.75 2.33
1.5.3 짝진 이표본 평균 - t.test(paired=TRUE)
- 귀무가설 : 모평균의 차이가 0이다.
- p-값이 작으므로 귀무가설 기각 즉, 두 방법에는 차이가 있다.
t.test(sleep$extra[sleep$group==1],sleep$extra[sleep$group==2],paired=TRUE)##
## Paired t-test
##
## data: sleep$extra[sleep$group == 1] and sleep$extra[sleep$group == 2]
## t = -4.0621, df = 9, p-value = 0.002833
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -2.4598858 -0.7001142
## sample estimates:
## mean of the differences
## -1.58
1.5.4 일표본 비율 1 - prop.test(x,n,p)
- 표본으로부터 모집단 비율 추정 및 가설 검정
- 정규 분포 근사 사용
- 신뢰구간에 0.5가 포함되어있고, p-value도 커서 귀무가설 채택 -> 찬성비율이 0.5라고 할 수 있다.
prop.test(42,100)##
## 1-sample proportions test with continuity correction
##
## data: 42 out of 100, null probability 0.5
## X-squared = 2.25, df = 1, p-value = 0.1336
## alternative hypothesis: true p is not equal to 0.5
## 95 percent confidence interval:
## 0.3233236 0.5228954
## sample estimates:
## p
## 0.42
prop.test(42,100,0.1)##
## 1-sample proportions test with continuity correction
##
## data: 42 out of 100, null probability 0.1
## X-squared = 110.25, df = 1, p-value < 2.2e-16
## alternative hypothesis: true p is not equal to 0.1
## 95 percent confidence interval:
## 0.3233236 0.5228954
## sample estimates:
## p
## 0.42
prop.test(42,100,0.5)##
## 1-sample proportions test with continuity correction
##
## data: 42 out of 100, null probability 0.5
## X-squared = 2.25, df = 1, p-value = 0.1336
## alternative hypothesis: true p is not equal to 0.5
## 95 percent confidence interval:
## 0.3233236 0.5228954
## sample estimates:
## p
## 0.42
1.5.5 일표본 비율 2 - binom.test(x,n,p)
- 표본으로부터 모집단 비율 추정 및 가설 검정
- 정규 분포 근사 사용하지 않고, 신뢰구간 직접 계산
binom.test(42,100)##
## Exact binomial test
##
## data: 42 and 100
## number of successes = 42, number of trials = 100, p-value = 0.1332
## alternative hypothesis: true probability of success is not equal to 0.5
## 95 percent confidence interval:
## 0.3219855 0.5228808
## sample estimates:
## probability of success
## 0.42
binom.test(42,100,0.5)##
## Exact binomial test
##
## data: 42 and 100
## number of successes = 42, number of trials = 100, p-value = 0.1332
## alternative hypothesis: true probability of success is not equal to 0.5
## 95 percent confidence interval:
## 0.3219855 0.5228808
## sample estimates:
## probability of success
## 0.42
1.5.6 이표본 비율 prop.test(c(x1,x2),c(n1,n2))
- 두 표본으로부터 모집단의 비율 추정 및 가설 검정에 사용
- 남성의 흡연율과 여성의 흡연율에 차이가 있을까?
- 남성 100명 중 45명, 여성 90명 중 55명 흡연한다면 두 비율에 차이가 있는가?
- 귀무가설 : 두 비율이 같다.
- p-value가 0.05보다 작으므로 귀무가설 기각 -> 두 비율에는 차이가 있다.
prop.test(c(45,55),c(100,90))##
## 2-sample test for equality of proportions with continuity
## correction
##
## data: c(45, 55) out of c(100, 90)
## X-squared = 4.3067, df = 1, p-value = 0.03796
## alternative hypothesis: two.sided
## 95 percent confidence interval:
## -0.31185005 -0.01037217
## sample estimates:
## prop 1 prop 2
## 0.4500000 0.6111111