如何在R中写出多项式回归?

bf1o4zei  于 2023-01-22  发布在  其他
关注(0)|答案(2)|浏览(165)

我正试图通过最小二乘法解出下面的公式:

数据:

dat <- structure(list(x = c(-3.06315207481384, -3.53966951370239, -0.786132514476776, 
2.85316586494446, 1.91431069374084, 0.238738358020782, -2.71875023841858, 
-4.51564025878906, -1.76425766944885, 36.0447616577148, 36.2231864929199, 
35.088565826416, 14.5400581359863, 11.0608978271484, 13.5029678344727, 
9.24409198760986, 9.52021884918213, 11.5315561294556, 9.64649677276611, 
-7.40514183044434, -5.99363708496094, -8.3693323135376, -6.55220413208008, 
16.6599369049072, 15.3896722793579, 14.0078706741333, -5.96392107009888, 
-6.05945920944214, -3.35051441192627, -4.49518489837646), y = c(-100.958953857422, 
-96.9687881469727, -99.6262283325195, -59.8251037597656, -61.4415702819824, 
-63.2435073852539, -99.5627517700195, -104.385040283203, -102.582824707031, 
-62.4078750610352, -60.3811149597168, -63.3329849243164, -87.6252822875977, 
-86.7032775878906, -89.111930847168, -48.4490776062012, -44.3725166320801, 
-44.0625152587891, -49.8619041442871, -99.2568817138672, -101.303993225098, 
-101.118591308594, -104.259971618652, -59.4827270507812, -61.4948043823242, 
-59.0172729492188, -114.305473327637, -113.537368774414, -115.28190612793, 
-115.021911621094), z = c(213.060043334961, 214.551696777344, 
214.275970458984, 186.614593505859, 189.293151855469, 185.896148681641, 
190.264511108398, 192.610168457031, 190.945388793945, 196.600631713867, 
197.723648071289, 197.753433227539, 179.977508544922, 180.026229858398, 
181.597991943359, 193.011459350586, 193.536712646484, 195.603088378906, 
194.673904418945, 212.374145507812, 212.009017944336, 210.889434814453, 
210.546630859375, 212.190979003906, 212.855712890625, 211.486724853516, 
203.421585083008, 205.626281738281, 206.354309082031, 202.603637695312
), activation_timing_ms = c(-8.95099963430584, -11.000197429963, 
-8.00000000000023, -29.3800032106763, -25.0499782952136, -29.9494419934181, 
-6.33357028643786, -7.84766776411197, -7, 4.41525051236772, 11, 
6.54840968699932, 8.55739954149135, -15.6434811032107, 15, -6.4484910424228, 
-7.04479069209742, -12.0007923079093, -8.9730733157603, -8.37403162067676, 
-9.07784949298457, -10.8114914004932, -6.9208993234206, -33.6635361890419, 
-27.2229083353382, -33.3012393050494, 1.82628266289839, 2, 0, 
1)), row.names = c(NA, 30L), class = "data.frame")

于是我写道:

lm <- lm(activation_timing_ms ~ x^3 + y^3 + z^3 + 
x^2*y + x^2*z + y^2*x + y^2*z + z^2*x + z^2*y + x*y*z + 
x^2 + y^2 + z^2 + x*y + x*z + y*z + x + y + z,
             data = dat)

这和我预期的不一样,因为它只返回了线性分量和交互作用。

coef(lm)

  (Intercept)             x             y             z           x:y 
-3.930285e+02  6.376611e+01 -4.131947e+00  1.754539e+00  5.855662e-01 
          x:z           y:z         x:y:z 
-3.100900e-01  1.895621e-02 -2.830448e-03

谢谢你的帮助!

wxclj1h5

wxclj1h51#

需要在lm中使用I()函数
从帮助文件中,您可以阅读:
在函数公式中,用于禁止将“+"、“-"、“*”和“^”等运算符解释为公式运算符,因此它们被用作算术运算符。它被terms. formula解释为一个符号。

lm(activation_timing_ms ~ I(x^3) + I(y^3) + I(z^3) + 
     I(x^2)*y + I(x^2)*z + I(y^2)*x + I(y^2)*z + I(z^2)*x + I(z^2)*y + x*y*z + 
     I(x^2) + I(y^2) + I(z^2) + x*y + x*z + y*z + x + y + z,
   data = dat)

Coefficients:
(Intercept)       I(x^3)       I(y^3)       I(z^3)       I(x^2)            y            z       I(y^2)            x  
 -5.514e+04   -1.263e-02    5.575e-04    5.903e-03    2.991e+00   -2.699e+01    8.109e+02    6.264e-01   -2.756e+02  
     I(z^2)     I(x^2):y     I(x^2):z     I(y^2):x     z:I(y^2)     x:I(z^2)     y:I(z^2)          y:x          z:x  
 -3.847e+00    1.291e-02   -8.344e-03   -6.022e-03   -2.320e-03   -2.948e-03   -3.120e-03   -2.483e+00    1.737e+00  
        y:z        y:z:x  
  8.326e-01    6.996e-03
uemypmqf

uemypmqf2#

你可以用

lm(activation_timing_ms~polym(x,y,z, degree=3, raw = TRUE), dat)

甚至

lm(activation_timing_ms~poly(cbind(x,y,z), 3, raw = TRUE), dat)

这导致:

summary(my_model)
Coefficients:
                                           Estimate Std. Error t value Pr(>|t|)  
(Intercept)                              -5.514e+04  4.439e+04  -1.242   0.2425  
poly(cbind(x, y, z), 3, raw = TRUE)1.0.0 -2.756e+02  3.265e+02  -0.844   0.4183  
poly(cbind(x, y, z), 3, raw = TRUE)2.0.0  2.991e+00  2.364e+00   1.266   0.2344  
poly(cbind(x, y, z), 3, raw = TRUE)3.0.0 -1.263e-02  1.430e-02  -0.883   0.3980  
poly(cbind(x, y, z), 3, raw = TRUE)0.1.0 -2.699e+01  2.217e+02  -0.122   0.9055  
poly(cbind(x, y, z), 3, raw = TRUE)1.1.0 -2.483e+00  1.246e+00  -1.993   0.0742 .
poly(cbind(x, y, z), 3, raw = TRUE)2.1.0  1.291e-02  1.639e-02   0.788   0.4492  
poly(cbind(x, y, z), 3, raw = TRUE)0.2.0  6.264e-01  4.625e-01   1.354   0.2055  
poly(cbind(x, y, z), 3, raw = TRUE)1.2.0 -6.022e-03  7.728e-03  -0.779   0.4539  
poly(cbind(x, y, z), 3, raw = TRUE)0.3.0  5.575e-04  1.575e-03   0.354   0.7307  
poly(cbind(x, y, z), 3, raw = TRUE)0.0.1  8.109e+02  6.611e+02   1.227   0.2481  
poly(cbind(x, y, z), 3, raw = TRUE)1.0.1  1.737e+00  2.891e+00   0.601   0.5614  
poly(cbind(x, y, z), 3, raw = TRUE)2.0.1 -8.344e-03  5.993e-03  -1.392   0.1940  
poly(cbind(x, y, z), 3, raw = TRUE)0.1.1  8.326e-01  2.294e+00   0.363   0.7241  
poly(cbind(x, y, z), 3, raw = TRUE)1.1.1  6.996e-03  8.067e-03   0.867   0.4061  
poly(cbind(x, y, z), 3, raw = TRUE)0.2.1 -2.320e-03  3.316e-03  -0.699   0.5002  
poly(cbind(x, y, z), 3, raw = TRUE)0.0.2 -3.847e+00  3.381e+00  -1.138   0.2818  
poly(cbind(x, y, z), 3, raw = TRUE)1.0.2 -2.948e-03  7.442e-03  -0.396   0.7003  
poly(cbind(x, y, z), 3, raw = TRUE)0.1.2 -3.120e-03  6.741e-03  -0.463   0.6534  
poly(cbind(x, y, z), 3, raw = TRUE)0.0.3  5.903e-03  5.973e-03   0.988   0.3463

相关问题