distribution - R - creating a linear model with fixed poly() coefficients -


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