앤디 필드의 유쾌한 R 통계학 챕터 6 - 상관분석

상관분석을 진행하고자 할 때 언제 어떤 방법을 사용할지 잘 판단해야 한다.

정규분포하는 연속형 변수 사이의 관계를 볼 때는 가장 일반적인 피어슨 상관계수를 사용한다.

정규분포가 아닌 변수 사이의 상관관계를 파악할 때는 스피어만 상관계수를 사용한다.

다만 중복값(동순위값)이 많은 경우에는 켄달의 타우를 사용한다.

 

연속변수와 이분변수 사이의 관계를 파악할 때 이연 상관계수를 사용하는데

이분변수의 원천이 연속형인 경우에는 이연, 아닌 경우에는 점이연이라고 하며

이연변수 확인을 위해서는 cor, cor.test 함수가 아닌 polyserial 함수를 사용한다.

(점이연은 그냥 cor함수를 사용한다.)

 

여러 변수들 중 내가 원하는 두 변수 사이의 관계만 확인하고 싶은 경우 편상관 분석을 시행한다.

이때 pcor 함수와 pcor.test 함수가 사용되는데 함수의 사용 방법이 cor, cor.test와 약간 다르니 주의해야 한다.

 

# 6장에서는 상관관계와 관련된 이론과 분석 방법을 배운다.

# 두 변수의 공분산을 계산하는 공식은 각 변수의 편차를 서로 곱한 것의 총합을 N-1로 나누는 것이다.
# 공분산이 양수라면 두 변수가 양의 상관관계를 갖고 음수라면 음의 상관관계를 갖는다고 할 수 있다.

# 여러 공분산을 비교하기 위해서는 각각의 척도를 표준화하여 사용한다.
# 편차를 표준편차로 나누면 편차가 표준화되어 몇 표준편차(시그마)만큼 떨어져 있는지 알 수 있다.

# 공분산을 표준화하기 위해서는 공분산을 두 변수의 표준편차의 곱으로 나눈다. 이를 (피어슨) 상관계수(r)라고 한다.
# 상관계수는 -1에서 1 사이의 값을 가지며 .1인 경우 작은 효과 .3인 경우 중간 효과, .5면 큰 효과로 판단한다.
# 상관계수가 유의한지를 판단하기 위해서 주로 자유도가 N-2인 t 통계량을 사용한다.

# 피어슨 상관계수를 사용하기 위해서는 1. 자료가 정규분포이고 2. 범주형 변수여야 한다.

# 산점도 그리기. 상관관계는 산점도를 통해 눈으로 확인할 수 있다.
adverts <- c(5,4,4,6,8)
packets <- c(8,9,10,13,15)
advertData <- data.frame(adverts, packets)

library(ggplot2)
a <- ggplot(advertData, aes(adverts, packets)) + geom_point()
a

# 필요한 패키지 불러오기(ggplot2는 이미 불러왔기 때문에 생략)
library(boot);library(ggm);library(Hmisc);library(polycor)

examData <- read.delim('datafiles/Exam Anxiety.dat', header = TRUE)

# 상관계수 도출을 위한 cor 함수의 사용
# use는 결측치 처리 방법으로 everything, all.obs, complete.obs, pairwise.complete.obs 중 선택하여 입력한다.
  # everything, all.obs는 데이터에 결측치가 존재하면 각각 NA, 오류메시지를 출력한다.
  # complete.obs는 결측치가 하나도 없는 행만 활용한다.
  # pairwise.complete.obs는 비교하는 두 열에 결측치가 없다면 해당 행 데이터를 활용한다.(나머지 열에 결측 있어도 활용)
# method는 상관계수의 종류로 pearson, spearman, kendall 중 선택하여 입력한다.

cor(examData$Exam,examData$Anxiety, use = 'complete.obs', method = 'pearson')

# 상관검정을 위한 cor.test 함수의 사용
  # alternative는 분석의 방향으로 two.sided, less, greater를 사용한다. 

cor.test(examData$Exam, examData$Anxiety, alternative = 'two.sided', method = 'pearson', conf.level = .95)

# 복습

examData2 <- examData[, c('Exam','Anxiety','Revise')]

cor(examData2, use = 'pairwise.complete.obs', method = 'pearson')

examMatrix <- as.matrix(examData2)

rcorr(examMatrix)

cor.test(examData$Exam, examData$Anxiety)


# 상관계수를 제곱한 것을 결정계수라고 하며 한 변수의 변동을 다른 변수가 얼마나 공유하는지를 표시한다.
# 결정계수가 0.153인 경우 y는 x의 변동을 15.3% 공유한다.

cor(examData2)^2

# 정규분포가 아닌 자료의 상관 분석을 위해서는 비모수적 통계인 스피어만 상관계수를 사용한다.

liarData <- read.delim('datafiles/The Biggest Liar.dat', header = TRUE)

cor(liarData$Position, liarData$Creativity, method = 'spearman')
cor.test(liarData$Position, liarData$Creativity, alternative = 'less', method = 'spearman')

# 동순위 점수가 많을 때에는 켄달의 타우를 사용한다.

cor(liarData$Position, liarData$Creativity, method = 'kendall')
cor.test(liarData$Position, liarData$Creativity, alternative = 'less', method = 'kendall')

# 부트스트랩을 통한 상관분석
# 객채 <- boot(자료, 함수, 반복)

bootTau <- function(liarData,i){cor(liarData$Position[i], liarData$Creativity[i], use = 'complete.obs', method = 'kendall')}

library(boot)
boot_kendall <- boot(liarData, bootTau, 2000)
boot_kendall
boot.ci(boot_kendall)


# 이연 상관과 점이연 상관
# 이연 상관은 연속 변수와 이분 변수지만 뒷배경에 어떤 연속체가 존재하는 변수 사이에서 사용한다.(ex. 시험 점수에 따른 합격/불합격)
# 점이연 상관은 연속 변수와 이산적인 이분 변수가 있을 때 사용한다.(ex. 죽음/생존)

catData <- read.csv('datafiles/pbcorr.csv', header = TRUE)

# 점이연 상관은 일반적인 상관 분석과 동일하게 진행할 수 있다.

cor.test(catData$time, catData$gender, method = 'pearson')

# 이연 상관은 polyserial 함수를 사용한다.

polyserial(catData$time, catData$gender)


# 편상관분석. 특정 변수의 효과를 통제하여 비교하고 싶은 변수 사이의 고유한 상관 관계를 분석한다.
# 편상관분석에는 pcor 함수와 pcor.test 함수를 사용한다.
# pcor(c('변수1', '변수2', '통제할 변수1', '통제할 변수 2'...), var(데이터프레임))
# pcor.test(pcor객체, 제어 변수의 수, 표본 크기)

pc <- pcor(c('Exam','Anxiety','Revise'), var(examData2));pc;pc^2
pcor.test(pc, 1, 103)


# 준편상관(부분 상관)
# 편상관분석에서는 여러 변수들 중 두 변수 사이의 관계에만 집중한다.
# 준편상관은 결과 변수의 변동을 다른 한 변수의 변동으로 설명하기 위해 통제하는 변수가 결과 변수에 미치는 영향을 제거한다. 편상관분석에서는 다른 모든 변수에 미치는 영향이 같이 제거되는 것과는 차이가 있다.

# 도출한 상관계수들 사이에 유의한 차이가 있는지 확인해보자
# 서로 다른 표본에 대한 상관계수의 차이는 z값을 통해 비교한다.
# 같은 표본에서 도출한 상관계수의 차이는 t값을 통해 비교한다.

# 복습 1
pfData <- read.delim('datafiles/Chamorro-Premuzic.dat', header = TRUE)
head(pfData)

cor(pfData[,c(3:7)], pfData[,c(8:12)], use = 'pairwise.complete.obs')

# 복습 2
essayData <- read.delim('datafiles/EssayMarks.dat', header = TRUE)
head(essayData)

ggplot(essayData, aes(hours, essay)) + geom_point()
ggplot(essayData, aes(hours, grade)) + geom_point()

cor.test(essayData$hours, essayData$essay)

# 복습 3
flickData <- read.delim('datafiles/ChickFlick.dat', header = TRUE)
head(flickData)

polyserial(flickData$arousal, flickData$gender)
polyserial(flickData$arousal, flickData$film)

# 복습 4
statData <- read.csv('datafiles/grades.csv', header = TRUE)
head(statData)

cor.test(statData$stats, statData$gcse, method = 'kendall')