parent
f57b81be21
commit
ee68dccbe9
@ -0,0 +1,85 @@
|
||||
rm(list=ls())
|
||||
set.seed(1123)
|
||||
|
||||
source('15_loocv.R')
|
||||
|
||||
multdiag <-
|
||||
function(X,d)
|
||||
{
|
||||
R <- matrix(NA, nrow=dim(X)[1], ncol=dim(X)[2])
|
||||
for (i in 1:dim(X)[2]) { R[,i]=X[,i]*d[i] }
|
||||
return(R)
|
||||
}
|
||||
|
||||
n <- 700
|
||||
p <- 55
|
||||
|
||||
sd <- 6 # standard deviation for zero-mean gaussian noise
|
||||
X <- matrix(runif(n*p),nrow=n,ncol=p)
|
||||
X <- scale(X)
|
||||
beta <- runif(p, min=-10, max=10)
|
||||
y <- X%*%beta + rnorm(n, mean=0, sd=sd)
|
||||
|
||||
lambdas <- 10^seq(-1,3,by=0.2)
|
||||
|
||||
var <-
|
||||
function(lambda)
|
||||
{
|
||||
d <- (Xs$d^2)/(Xs$d^2 + lambda)^2
|
||||
var <- multdiag(Xs$v,d)
|
||||
var <- sd^2 * tcrossprod(var,Xs$v)
|
||||
}
|
||||
|
||||
bias <-
|
||||
function(lambda)
|
||||
{
|
||||
d <- lambda/(Xs$d^2+lambda)
|
||||
bias <- multdiag(Xs$v,d)
|
||||
bias <- bias %*% crossprod(Xs$v,beta)
|
||||
}
|
||||
|
||||
epeVar <-
|
||||
function(lambda)
|
||||
{
|
||||
var <- var(lambda)
|
||||
return( mean(rowSums(X*(X%*%var))) )
|
||||
}
|
||||
|
||||
epeBias <-
|
||||
function(lambda)
|
||||
{
|
||||
bias <- bias(lambda)
|
||||
return( mean((X%*%bias)^2) )
|
||||
}
|
||||
|
||||
epe <-
|
||||
function(lambda)
|
||||
{
|
||||
return( epeVar(lambda) + epeBias(lambda) + sd^2 )
|
||||
}
|
||||
|
||||
#rm <- ridge(X, y, lambdas)
|
||||
#Xs <- svd(X)
|
||||
#epes <- sapply(lambdas, epe)
|
||||
|
||||
X.init <- X
|
||||
beta.init <- beta
|
||||
|
||||
ps <- seq(1,p)
|
||||
lps <- length(ps)
|
||||
epes <- numeric(lps)
|
||||
biass <- numeric(lps)
|
||||
vars <- numeric(lps)
|
||||
maes <- numeric(lps)
|
||||
|
||||
for (k in ps)
|
||||
{
|
||||
X <- X.init[,1:k]
|
||||
beta <- beta.init[1:k]
|
||||
rm <- ridge(X, y, lambdas)
|
||||
Xs <- svd(X)
|
||||
epes[k] <- epe(rm$lambda)
|
||||
biass[k] <- epeBias(rm$lambda)
|
||||
vars[k] <- epeVar(rm$lambda)
|
||||
maes[k] <- rm$mae
|
||||
}
|
@ -0,0 +1,4 @@
|
||||
#!/bin/bash
|
||||
|
||||
arg="library(bookdown);bookdown::preview_chapter(\"$1\")"
|
||||
Rscript -e $arg
|
Loading…
Reference in new issue