гладкий.spline (): установленная модель не соответствует заданной пользователем степени свободы

Вот код, который я запустил 

fun <- function(x) {1 + 3*sin(4*pi*x-pi)}
set.seed(1)
num.samples <- 1000
x <- runif(num.samples)
y <- fun(x) + rnorm(num.samples) * 1.5
fit <- smooth.spline(x, y, all.knots=TRUE, df=3)

Несмотря на df=3, Когда я проверил подходящую модель, выход был 

Call:
smooth.spline(x = x, y = y, df = 3, all.knots = TRUE)
Smoothing Parameter  spar= 1.499954  lambda= 0.002508571 (26 iterations)
Equivalent Degrees of Freedom (Df): 9.86422

Не мог бы кто-нибудь помочь? Спасибо!

1 5
regression

1 ответ:

. отметим, что из Р-3.4.0 (2017-04-21), smooth.spline можете принять прямое указание λ недавно добавил аргумент lambda. Но он все равно будет преобразован во внутренний spar во время оценки. Таким образом, следующий ответ не влияет.

Параметр сглаживания λ / spar лежит в центре управления гладкостью

Гладкость контролируется параметром сглаживания λ.smooth.spline() использует внутренний параметр сглаживания spar, а не λ: 
spar = s0 + 0.0601 * log(λ)
Такое логарифмическое преобразование необходимо для того, чтобы сделать неограниченную минимизацию, как GCV/CV. Пользователь может указать spar, чтобы косвенно указать λ. Когда spar растет линейно, λ будет расти экспоненциально. Таким образом, редко возникает необходимость в использовании большого значения spar.

Степень свободы df, также определяется в терминах λ:

ЭДФ

Где X - матрица модели с базисом B-сплайна и S - Штраф матрица.

Вы можете проверить их отношения с вашим набором данных: 

spar <- seq(1, 2.5, by = 0.1)
a <- sapply(spar, function (spar_i) unlist(smooth.spline(x, y, all.knots=TRUE, spar = spar_i)[c("df","lambda")]))

Давайте набросаем df ~ spar, λ ~ spar и log(λ) ~ spar: 

par(mfrow = c(1,3))
plot(spar, a[1, ], type = "b", main = "df ~ spar",
     xlab = "spar", ylab = "df")
plot(spar, a[2, ], type = "b", main = "lambda ~ spar",
     xlab = "spar", ylab = "lambda")
plot(spar, log(a[2,]), type = "b", main = "log(lambda) ~ spar",
     xlab = "spar", ylab = "log(lambda)")

сюжет

Обратите внимание на радикальный рост λ с spar, линейную связь между log(λ) и spar и относительно гладкую связь между df и spar.

smooth.spline() подходящие итерации для spar

Если мы вручную зададим значение spar, как то, что мы сделали в sapply() нет подходящих итераций для выбора spar; в противном случае smooth.spline() необходимо перебрать ряд значений spar. Если мы

  • укажите cv = TRUE / FALSE, подгонка итераций направлена на минимизацию CV / GCV Баллов;
  • укажите df = mydf, подгонка итераций направлена на минимизацию (df(spar) - mydf) ^ 2.

За минимизацией GCV легко следить. Мы не заботимся о счете GCV, но заботимся о соответствующем spar. Напротив, при минимизации (df(spar) - mydf)^2 мы часто заботимся о df значение в конце итерации, а не spar! Но принимая во внимание, что это проблема минимизации, мы никогда не гарантируем, что конечное df соответствует нашему целевому значению mydf.

Почему вы ставите df = 3, но получаете df = 9.864?

Окончание итерации может означать либо достижение минимума, либо достижение границы поиска, либо достижение максимального числа итераций.

Мы далеки от максимального предела итераций (по умолчанию 500); тем не менее мы не попали в минимальный. Ну, мы могли бы достичь границы.

Не сосредотачивайтесь на df, думайте о spar. 
smooth.spline(x, y, all.knots=TRUE, df=3)$spar   # 1.4999

В соответствии с ?smooth.spline, по умолчанию, smooth.spline() ищет spar между [-1.5, 1.5]. То есть, когда вы ставите df = 3, минимизация заканчивается на границе поиска, а не на df = 3.

Еще раз взгляните на наш график отношения между df и spar. Из рисунка видно, что нам нужно некоторое значение spar около 2, чтобы привести к [50]}.

Воспользуемся аргументом control.spar: 

fit <- smooth.spline(x, y, all.knots=TRUE, df=3, control.spar = list(high = 2.5))
# Smoothing Parameter  spar= 1.859066  lambda= 0.9855336 (14 iterations)
# Equivalent Degrees of Freedom (Df): 3.000305
Теперь вы видите, что вы заканчиваете с df = 3. И нам нужен spar = 1.86.

Лучшее предложение: не используйте all.knots = TRUE

Посмотрите, у вас есть 1000 данных. С помощью all.knots = TRUE вы будете использовать 1000 параметров. Желание закончить с df = 3 подразумевает, что 997 из 1000 параметров подавляются. Представьте, какой большой a λ следовательно spar вам нужен!

Попробуйте вместо этого использовать сплайн штрафной регрессии. Подавление 200 параметры до 3 определенно намного проще: 
fit <- smooth.spline(x, y, nknots = 200, df=3)  ## using 200 knots
# Smoothing Parameter  spar= 1.317883  lambda= 0.9853648 (16 iterations)
# Equivalent Degrees of Freedom (Df): 3.000386
Теперь вы заканчиваете с df = 3 Без spar контроля.
5