in r, 1 can build lm()
or glm()
object fixed coefficients, using offset
parameter in formula.
x=seq(1,100) y=x^2+3*x+7 # forcing fit polynomial: 2x^2 + 4x + 8 fixed_model = lm(y ~ 0 + offset(8 + 4*x + 2*i(x^2) ))
is possible same thing using poly()
? tried code below doesn't seem work.
fixed_model_w_poly <- lm(y ~ 0 + offset(poly(x, order=2, raw = true, coefs= c(8, 4, 2))))
error : number of offsets 200, should equal 100 (number of observations)
i want use poly()
convenient interface run iterations high number of fixed coefficients or order values, rather having manually code: offset(8 + 4*x + 2*i(x^2) )
each order/coefficient combination.
p.s: further not essential information: go inside mcmc routine. example usage generate (and compare) model_current
model_next
in below code:
library(mass) coeffs_current <- c(8, 4, 2) model_current <- lm(y ~ 0 + offset(poly(x, order=2, raw = true, coefs= coeffs_current ))) cov <- diag(rep(1,3)) coeffs_next <- mvrnorm(1, mu = as.numeric(coeffs_current ), sigma = cov ) model_next <- lm(y ~ 0 + offset(poly(x, order=2, raw = true, coeffs_next ))
this demonstrates suggested. (not use poly.)
library(mass) # coeffs_current <- c(8, 4, 2) name change compactness. cc <- c(8, 4, 2) form <- as.formula(bquote(y~x+offset(.(cc[1])+x*.(cc[2])+.(cc[3])*i(x^2) ))) model_current <- lm(form, data=dat))
i have no idea intend next code. looks want based on inputs prior function, doesn't want based on results.
cov <- diag(rep(1,3)) coeffs_next <- mvrnorm(1, mu = as.numeric(cc ), sigma = cov )
the code works (at least intended) simple test case. bquote function substitutes values expressions (well calls) , as.formula function evaluates argument , dresses result proper formula
-object.
dat <- data.frame(x=rnorm(20), y=rnorm(20) ) cc <- c(8, 4, 2) form <- as.formula( bquote(y~x+offset(.(cc[1])+x*.(cc[2])+.(cc[3])*i(x^2) ))) model_current <- lm(form, data=dat) #-------- > model_current call: lm(formula = form, data = dat) coefficients: (intercept) x -9.372 -5.326 # bizarre results due offset. #-------- form #y ~ x + offset(8 + x * 4 + 2 * i(x^2))
Comments
Post a Comment