통계 이야기

data mining, polynomial regression + step functions + Natural cubic spline + smoothing spline + local regression + GAM in R

창이 2021. 12. 18.
728x90
반응형

polynimial regression (다항식 회귀) 

- y_i=β_0+β_1 x_i+β_2 x_i^2_3 x_i^3+〖…+β〗_d x_i^d+ϵ_i   [ϵ_i is the error term]
- Generally, d is not greater than 3 or 4 (더 커지면 너무 극심하게 비선형 곡선이 됨)

step functions (계단 함수)

- X의 범위를 여러 개의 bin으로 분할하여 각 bin에 다른 상수를 적합
- Continuous variableordered categorical variable로 변환

 

regression splines (회귀 스플라인)

- Piecewise polynomials regression with a single knot (단일 매듭 조각별 다항식 회귀)
1.Piecewise cubic
- x_i<c인 관측치로 첫번째 cubic을 적합 시키고 x_ic 관측치로 두번째 cubic을 적합
- Parameter개수가 8개 이므로 적합 시에 df=8이 사용됨
2. Continuous piecewise cubic
- Piecewise cubic에서 한 개의 constraint 추가: continuity (lim┬(x_i→c)⁡〖f(x_i)〗=f(c))
- 한 개의 constraint로 자유도 df=8-1=7이 사용됨. 연속적이지만 자연스럽지는 않음
3.Cubic spline
- Continuous piecewise cubic에서 두 개의 constraint 추가: continuity of the first and second derivatives (lim┬(x_i→c)⁡
- 두개의 constraint 추가로 자유도 df=8-1-2=5가 사용됨.
4. Linear spline: linear function
두개로 적합하고 연속성 조건을 만족시키면 됨

Natural cubic spline (자연 삼차 스플라인)
- spline의 양쪽 경계에서 선형이라는 추가적인 constraint가 있어서 더 안정적인 추정치를 제공
Number and Locations of the Knots
 
- Number of the knots
 
- 10-fold CV를 통해서 최소의 test MSE를 주는 K값을 선택 (아래 그림에서 K=3,4면 충분)
 
Location of the knots
 
- knot의 개수가 결정되면 uniform quantile에 대응하는 수의 knot을 위치 시킴
 
- 3개의 knot인 경우 age25th, 50th, 75th percentileknot의 위치로 사용

smoothing splines

- training observations: x_1,…, x_n
- tuning parameter λ controls the roughness of the smoothing spline
- effective degree of freedom (df_λ): λ increases from 0 to df_λ decreases n to 2 (직선)
- leave-one-out CV (LOOCV)RSS_cv가 최소가 되게 하는 λ값을 선택

local regression(국소 회귀)

- compute the fit at a target point x_0 using only the regression nearby training observations
- Weight function K: 다양한 방법으로 정의 가능하며 linear, quadratic으로 적합 가능
- Span s = k/n
- smoothing splinetuning parameter처럼 non-linear fitflexibility를 조절함. s가 작을 수록 적합은 더 local하고 꾸불꾸불해지고 s가 아주 크면 training data 전체를 사용하여 global fit을 하게 
- K-fold CV, LOOCV 등의 방법으로 s를 선택할 수 있음
- 그림에서 s=0.2, 0.7로 선택하여 적합. s=0.7 때가 s=0.2일 때 보다 smooth

 

GAM for regression problems

- explore the problem of flexibly predicting Y on the basis of several predictors, X_1, . . . , X_p
- an extension of multiple linear regression
- extend a standard linear model by allowing non-linear functions of each of the variables, while maintaining additivity
- GAM의 장점은 linear model이나 non-linear model들을 additive model을 적합하기 위한 building block으로 사용할 수 있다는 점임
- GAM의 단점은 modeladditive해야 한다는 점임. 하지만 Linear regression처럼 interaction termGAM에 추가할 수 있음

# fit polynomial regression
fit=lm(wage~poly(age,4),data=Wage) # orthogonal polynomial
coef(summary(fit))
fit2=lm(wage~poly(age,4,raw=T),data=Wage)
coef(summary(fit2))
fit2a=lm(wage~age+I(age^2)+I(age^3)+I(age^4),data=Wage)
coef(fit2a)
fit2b=lm(wage~cbind(age,age^2,age^3,age^4),data=Wage)
coef(fit2b)
  # make plot for degree-4 polynomial regression
  
  agelims=range(age)
  age.grid=seq(from=agelims[1],to=agelims[2])
  preds=predict(fit,newdata=list(age=age.grid),se=TRUE)
  se.bands=cbind(preds$fit+2*preds$se.fit,preds$fit-2*preds$se.fit)
  
  par(mfrow=c(1,1),mar=c(4.5,4.5,1,1),oma=c(0,0,4,0))
  plot(age,wage,xlim=agelims,cex=.5,col="darkgrey")
  title("Degree-4 Polynomial",outer=T)
  lines(age.grid,preds$fit,lwd=2,col="blue")
  matlines(age.grid,se.bands,lwd=1,col="blue",lty=3)
  preds2=predict(fit2,newdata=list(age=age.grid),se=TRUE)
  max(abs(preds$fit-preds2$fit))

  # Select the degree of the polynomial
  fit.1=lm(wage~age,data=Wage)
  fit.2=lm(wage~poly(age,2),data=Wage)
  fit.3=lm(wage~poly(age,3),data=Wage)
  fit.4=lm(wage~poly(age,4),data=Wage)
  fit.5=lm(wage~poly(age,5),data=Wage)
  anova(fit.1,fit.2,fit.3,fit.4,fit.5);    coef(summary(fit.5))
  fit=glm(I(wage>250)~poly(age,4),data=Wage,family=binomial)
  preds=predict(fit,newdata=list(age=age.grid),se=T)
  pfit=exp(preds$fit)/(1+exp(preds$fit))
  se.bands.logit = cbind(preds$fit+2*preds$se.fit, preds$fit-2*preds$se.fit)
  se.bands = exp(se.bands.logit)/(1+exp(se.bands.logit))
  preds=predict(fit,newdata=list(age=age.grid),type="response",se=T)
  
  plot(age,I(wage>250),xlim=agelims,type="n",ylim=c(0,.2))
  points(jitter(age), I((wage>250)/5),cex=.5,pch="|",
         col="darkgrey")
  lines(age.grid,pfit,lwd=2, col="blue")
  matlines(age.grid,se.bands,lwd=1,col="blue",lty=3)

# step function
table(cut(age,4))
fit=lm(wage~
cut(age,4),data=Wage)
coef(summary(fit))
## (2) Splines

# regression splines
library(splines)
fit=lm(wage~bs(age,knots=c(25,40,60)),data=Wage) # B-Spline Basis
pred=predict(fit,newdata=list(age=age.grid),se=T)
plot(age,wage,col="gray")
lines(age.grid,pred$fit,lwd=2)
lines(age.grid,pred$fit+2*pred$se,lty="dashed")
lines(age.grid,pred$fit-2*pred$se,lty="dashed")
dim(bs(age,knots=c(25,40,60)))
dim(bs(age,df=6))
attr(bs(age,df=6),"knots")

# natural splines
fit2=lm(wage~ns(age,df=4),data=Wage) # natural splines
pred2=predict(fit2,newdata=list(age=age.grid),se=T)
lines(age.grid, pred2$fit,col="red",lwd=2)

plot(age,wage,xlim=agelims,cex=.5,col="darkgrey")
title("Smoothing Spline")
fit=smooth.spline(age,wage,df=16)
fit2=smooth.spline(age,wage,cv=TRUE)
fit2$df
lines(fit,col="red",lwd=2)
lines(fit2,col="blue",lwd=2)
legend("topright",legend=c("16 DF","6.8 DF"),
col=c("red","blue"),lty=1,lwd=2,cex=.8)

# local regression
plot(age,wage,xlim=agelims,cex=.5,col="darkgrey")
title("Local Regression")
fit=loess(wage~age,span=.2,data=Wage)
fit2=loess(wage~age,span=.5,data=Wage)
lines(age.grid,predict(fit,data.frame(age=age.grid)),
col="red",lwd=2)
lines(age.grid,predict(fit2,data.frame(age=age.grid)),
col="blue",lwd=2)
legend("topright",legend=c("Span=0.2","Span=0.5"),
col=c("red","blue"),lty=1,lwd=2,cex=.8)
## (3) GAMs

install.packages("gam")
library(gam)

# Fit smoothing splines
gam1=lm(wage~ns(year,4)+ns(age,5)+education,data=Wage)
gam.m3=gam(wage~s(year,4)+s(age,5)+education,data=Wage)
par(mfrow=c(1,3))
plot(gam.m3, se=TRUE,col="blue")

# Compare three models
gam.m1=gam(wage~s(age,5)+education,data=Wage)
gam.m2=gam(wage~year+s(age,5)+education,data=Wage)
anova(gam.m1,gam.m2,gam.m3,test="F")
summary(gam.m3)

 

728x90
반응형

댓글

추천 글