multivariate adaptive regression splines (MARS)

MARSは多変量非線形モデルです。

ライブラリは"mda"で、Hastie先生とTibshirani先生が作ってます（この先生達いろいろ作ってるな、、）。

サンプルデータはtreeというデータで結果変数はVolumeという変数。

データはこんな感じ。

library(mda)	#MARS: Trevor Hastie & Robert Tibshirani
library(gclus)

data(trees)
state.cor <- cor(trees)
cpairs(trees, panel.colors=dmat.color(state.cor), order.single(state.cor),pch=".",gap=.4)

VolumeはGirthと相関が高く、単純な直線関係ではなくてちょっと折れ曲がってる感じがする。

そこで2種類のMARSを行ってみます。

下記の関数はexampleのものを修正したものです。

#---calculate cut point
showcuts <- function(obj){
	tmp           <- obj$cuts[obj$sel, ]
	dimnames(tmp) <- list(NULL)
	tmp
}

#---draw fit curv
PlotMars <- function(XData, YVec, PredFit, ShowCuts = "No"){
	for(i in 1:ncol(XData)) {
		Xp       <- matrix(sapply(XData, mean), nrow(XData), ncol(XData), byrow=TRUE)
		xr       <- sapply(XData, range)
		Xp1      <- Xp
		Xp1[, i] <- seq(xr[1, i], xr[2, i], len = nrow(XData))
		Xf       <- predict(PredFit, Xp1)
		plot(XData[, i], YVec, xlab = names(XData)[i], ylab = "")
		lines(Xp1[ ,i], Xf, col = "blue")
		if(ShowCuts == "Yes"){
			abline(v = max(showcuts(PredFit)[, i]), lty = 2)
			text(max(showcuts(PredFit)[, i]), max(axTicks(2)) * 0.9, max(showcuts(PredFit)[, i]), adj = 0, col = "blue", cex = 2)
		}
	}
	RSS <- sum( (YVec - mean(YVec) )^2)
	Res <- as.numeric(PredFit$residuals)%*%as.numeric(PredFit$residuals) 
	Rsq <- (RSS - Res)/RSS
	text(min(axTicks(1)), max(axTicks(2)), paste("R-sq = ", round(Rsq, digits = 2), sep = ""), adj = 0)
}


fit1 <- bruto(trees[,-3], trees[3])
fit2 <- mars(trees[,-3], trees[3])

PlotMars(trees[,-3], unlist(trees[3]), fit1)
PlotMars(trees[,-3], unlist(trees[3]), fit2, ShowCuts = "Yes")

2つの図はmarsを行った結果。

Girth=12の点で回帰直線が折れ曲がっています。

bruto出の結果は折れ線じゃなくて曲線になっていて、こっちの方がfittingは良いみたいです。

ちなみに線形回帰だとこうなります。

summary(lm(Volume ~., data = trees))

R2乗は0.944なのでmarsの方が良さそう。

でもこういう非線形回帰はoverfittingが起こりやすいので、クロスバリデーションした方が良いですね。