#title 변수간의 상관정도에 대한 분석 [[TableOfContents]] ==== 상관분석(Correlation Analysis) ==== * 두 확률 변수의 상관 정도를 분석 (원인이 아님) * 가계소득이 높아지면 자녀의 성적도 높아지는가? * 공분산 값 자체의 의미보다는 부호가 중요한 의미를 가짐 (공분산은 방향성만 파악) * 피어슨 상관계수는 방향 및 정도에 대한 2가 의미 (두 변수가 모두 정규분포이어야..) * 스피어만, 켄달 상관계수는 분포를 모를때 사용 (켄달 상관계수가 모집단의 추정치이므로 더 선호된다고 함) * 상관분석 * 산점도(scatter plot) * 상관계수(correlation coefficient) ==== 예제데이터 ==== 허리둘레(girth)와 몸무게(weigth) {{{ tmp <- textConnection( "girth weight 35.00 3.45 32.00 3.20 30.00 3.00 31.50 3.20 32.70 3.30 30.00 3.20 36.00 3.85 30.50 3.15 34.70 3.65 30.50 3.40 33.00 3.50 35.00 4.00 31.80 3.10 38.00 4.20 33.00 3.45") x <- read.table(tmp, header=TRUE) close.connection(tmp) plot(x) }}} 산점도 attachment:변수간의상관정도에대한분석/1.png ==== 공분산(Covariance) ==== {{{ > cov(x) girth weight girth 5.7183810 0.7461667 weight 0.7461667 0.1210238 }}} 공분산이 양수이므로 양의 상관관계를 가짐 ==== 피어슨 상관계수 ==== {{{ > cor.test(girth, weight, data=x, method="pearson") #default: pearson Pearson's product-moment correlation data: girth and weight t = 7.3142, df = 13, p-value = 5.881e-06 alternative hypothesis: true correlation is not equal to 0 95 percent confidence interval: 0.7116680 0.9655589 sample estimates: cor 0.896941 }}} * 피어슨 상관계수가 0.896941로 강한 양의 상관관계가 있다고 볼 수 있음 * [http://ko.wikipedia.org/wiki/%EC%83%81%EA%B4%80%EB%B6%84%EC%84%9D#.ED.94.BC.EC.96.B4.EC.8A.A8_.EC.83.81.EA.B4.80.EA.B3.84.EC.88.98_.28Pearson_correlation_coefficient.29 해석] * 1.0 ~ 0.7: 강 * 0.7 ~ 0.3: 중 * 0.3 ~ 0.1: 약 * 0.1 ~ 0.0: 무시 ==== 스피어만 상관계수 ==== {{{ > cor.test(girth, weight, data=x, method = "spearman") Spearman's rank correlation rho data: girth and weight S = 70.0628, p-value = 1.963e-05 alternative hypothesis: true rho is not equal to 0 sample estimates: rho 0.8748878 Warning message: In cor.test.default(girth, weight, data = x, method = "spearman") : Cannot compute exact p-values with ties }}} * 순위를 가지고 하는 것이므로 같은 값이 있어 warning이 보임. * 상관계수가 0.8748878로 강한 양의 상관관계 참고: [attachment:변수간의상관정도에대한분석/spearman_corr.xlsx 엑셀 계산] ==== 켄달 상관계수 ==== {{{ > cor.test(girth, weight, data=x, method = "kendall") Kendall's rank correlation tau data: girth and weight z = 3.7508, p-value = 0.0001762 alternative hypothesis: true tau is not equal to 0 sample estimates: tau 0.7425743 Warning message: In cor.test.default(girth, weight, data = x, method = "kendall") : ties 때문에 정확한 p값을 계산할 수가 없습니다 }}} * 상관계수가 0.7425743로 역시 강한 상관관계 임. ==== data frame에서 상관계수 다루기 ==== {{{ install.packages("Hmisc") library("Hmisc") rcorr(as.matrix(x), type="pearson") rcorr(as.matrix(x), type="spearman") }}} ==== 시각화1 ==== {{{ #일반적인 산점도 pairs(x1) install.packages("GGally") library("GGally") #정보가 많은 산점도 ggpairs(x1) }}} ==== 시각화2 ==== 출처: R Graphics Cookbook {{{ panel.cor <- function(x, y, digits=2, prefix="", cex.cor, ...) { usr <- par("usr") on.exit(par(usr)) par(usr = c(0, 1, 0, 1)) r <- abs(cor(x, y, use="complete.obs")) txt <- format(c(r, 0.123456789), digits=digits)[1] txt <- paste(prefix, txt, sep="") if(missing(cex.cor)) cex.cor <- 0.8/strwidth(txt) text(0.5, 0.5, txt, cex = cex.cor * (1 + r) / 2) } panel.hist <- function(x, ...) { usr <- par("usr") on.exit(par(usr)) par(usr = c(usr[1:2], 0, 1.5) ) opar=par(ps=30) #font-size h <- hist(x, plot = FALSE) breaks <- h$breaks nB <- length(breaks) y <- h$counts y <- y/max(y) rect(breaks[-nB], 0, breaks[-1], y, col="white", ...) } panel.lm <- function (x, y, col = par("col"), bg = NA, pch = par("pch"), cex = 1, col.smooth = "black", ...) { points(x, y, pch = pch, col = col, bg = bg, cex = cex) abline(stats::lm(y ~ x), col = col.smooth, ...) #로버스한거 #abline(MASS::rlm(y ~ x), col = col.smooth, ...) } pairs(c2009[,2:5], pch=".", upper.panel = panel.cor, diag.panel = panel.hist, lower.panel = panel.lm) }}} ==== sql로 상관계수 구하기 ==== {{{ select (count(*) * sum(x * y) - sum(x) * sum(y)) / (sqrt(count(*) * sum(x * x) - sum(x) * sum(x)) * sqrt(count(*) * sum(y * y) - sum(y) * sum(y))) r from (select x,y from data) t }}} ==== 로버스트 상관 분석 ==== {{{ library("robust") covRob(stackloss) #Robust Estimate of Covariance covRob(stackloss, corr = TRUE) #Robust Estimate of Correlation }}} 결과 {{{ > covRob(stackloss, estim = "mcd") Call: covRob(data = stackloss, estim = "mcd") Robust Estimate of Covariance: Air.Flow Water.Temp Acid.Conc. stack.loss Air.Flow 77.32 24.02 54.91 67.30 Water.Temp 24.02 18.28 18.40 25.16 Acid.Conc. 54.91 18.40 99.03 48.99 stack.loss 67.30 25.16 48.99 60.93 Robust Estimate of Location: Air.Flow Water.Temp Acid.Conc. stack.loss 56.15 20.23 85.38 13.15 > covRob(stackloss, corr = TRUE) Call: covRob(data = stackloss, corr = TRUE) Robust Estimate of Correlation: Air.Flow Water.Temp Acid.Conc. stack.loss Air.Flow 1.0000 0.6677 0.6174 0.9514 Water.Temp 0.6677 1.0000 0.4960 0.7868 Acid.Conc. 0.6174 0.4960 1.0000 0.5389 stack.loss 0.9514 0.7868 0.5389 1.0000 Robust Estimate of Location: Air.Flow Water.Temp Acid.Conc. stack.loss 56.92 20.43 86.29 13.73 }}} ==== 편상관분석 ==== 두 변수간 상관관계에 영향을 줄 수 있는 제 3의 변수가 존재할 수도 있으므로 이를 통제하고 상관분석 수행 {{{ library(ggm) height <- c(160, 155, 170, 160, 180, 177) weight <- c(55, 50, 56, 70,80,73) gender <- c(0,0,0,1,1,1) mydata <- data.frame(height, weight, gender) pcor(c("height", "weight", "gender"), var(mydata)) # partial corr between [height] and [weight] controlling for [gender] cor(height, weight) }}} 결과 {{{ > pcor(c("height", "weight", "gender"), var(mydata)) [1] 0.8267963 > cor(height, weight) [1] 0.7598541 > }}} --http://wolfpack.hnu.ac.kr/2015_Fall/MDA/MDA_correlation.pdf {{{ ds=read.csv("SMSA_USA.csv") #데이터 읽기 attach(ds) #데이터 사용 선언; names(ds) #변수이름 리스트 library(car) #행렬 산점도 패키지 설치 scatterplot.matrix(~d_index+edu+non_white+wc_ratio+income) library(Hmisc) #상관계수 행렬 계산 패키지 설치 rcorr(as.matrix(ds[6:13]),type=“spearman") #부분결정계수 summary(lm(d_index~income)) #사망률(Y)=소득(X) 회귀분석 fit=lm(d_index~edu) # 회귀분석 Y=Z(통제변수, 교육수준) fit2=lm(income~edu) # 회귀분석 X=Z plot(fit$residuals, fit2$residuals) #산점도 (X, Y) : Z 통제 summary(lm(fit$residuals~fit2$residuals)) #통제 후 Y=X 회귀분석 cor.test(fit$residuals,fit2$residuals) #부분상관계수 }}} ==== 참고자료 ==== * [정준상관분석] --> 독립변수군과 종속변수군의 상관관계 ---- 오 좋네요 -- sion 2018-09-12 20:55:04