Tỷ lệ đa thức xấp xỉ

Question

Tôi cố gắng để phù hợp với một đa thức cho bộ dữ liệu của tôi, trông như thế (full bộ dữ liệu là ở phần cuối của bài viết):

Lý thuyết dự đoán rằng việc xây dựng các đường cong :

mà trông như thế này (cho x nằm giữa 0 và 1):

Khi tôi cố gắng để tạo ra một mô hình tuyến tính trong R bằng cách thực hiện:

mod <- lm(y ~ poly(x, 2, raw=TRUE)/poly(x, 2)) 

tôi nhận được đường cong sau:

Đó là khác nhiều so với những gì tôi mong đợi. Bạn có bất kỳ ý tưởng làm thế nào để phù hợp với một đường cong mới từ dữ liệu này để nó sẽ tương tự như một, mà lý thuyết dự đoán? Ngoài ra, nó chỉ nên có một mức tối thiểu.

Full bộ dữ liệu:

---

Vector các giá trị x:

x <- c(0.02, 0.03, 0.04, 0.05, 0.06, 0.07, 0.08, 0.09, 0.10, 0.11, 0.12, 
0.13, 0.14, 0.15, 0.16, 0.17, 0.18, 0.19, 0.20, 0.21, 0.22, 0.23, 0.24, 0.25, 
0.26, 0.27, 0.28, 0.29, 0.30, 0.31, 0.32, 0.33, 0.34, 0.35, 0.36, 0.37, 0.38, 
0.39, 0.40, 0.41, 0.42, 0.43, 0.44, 0.45, 0.46, 0.47, 0.48, 0.49, 0.50, 0.51, 
0.52, 0.53, 0.54, 0.55, 0.56, 0.57, 0.58, 0.59, 0.60, 0.61, 0.62, 0.63, 0.64, 
0.65, 0.66, 0.67, 0.68, 0.69, 0.70, 0.71, 0.72, 0.73, 0.74, 0.75, 0.76, 0.77, 
0.78, 0.79, 0.80, 0.81, 0.82, 0.83, 0.84, 0.85, 0.86, 0.87, 0.88, 0.89, 0.90, 
0.91, 0.92, 0.93, 0.94, 0.95) 

Vector các giá trị y:

y <- c(4.104, 4.444, 4.432, 4.334, 4.285, 4.058, 3.901, 4.382, 
    4.258, 4.158, 3.688, 3.826, 3.724, 3.867, 3.811, 3.550, 3.736, 3.591, 
    3.566, 3.566, 3.518, 3.581, 3.505, 3.454, 3.529, 3.444, 3.501, 3.493, 
    3.362, 3.504, 3.365, 3.348, 3.371, 3.389, 3.506, 3.310, 3.578, 3.497, 
    3.302, 3.530, 3.593, 3.630, 3.420, 3.467, 3.656, 3.644, 3.715, 3.698, 
    3.807, 3.836, 3.826, 4.017, 3.942, 4.208, 3.959, 3.856, 4.157, 4.312, 
    4.349, 4.286, 4.483, 4.599, 4.395, 4.811, 4.887, 4.885, 5.286, 5.422, 
    5.527, 5.467, 5.749, 5.980, 6.242, 6.314, 6.587, 6.790, 7.183, 7.450, 
    7.487, 8.566, 7.946, 9.078, 9.308, 10.267, 10.738, 11.922, 12.178, 13.243, 
    15.627, 16.308, 19.246, 22.022, 25.223, 29.752) 

Answer 1

Sử dụng nls để phù hợp với mô hình phi tuyến. Lưu ý rằng công thức mô hình không được định nghĩa duy nhất như được hiển thị trong câu hỏi vì nếu chúng ta nhân tất cả các hệ số với bất kỳ số nào, kết quả sẽ vẫn đưa ra các dự đoán tương tự. Để tránh điều này, chúng ta cần sửa một hệ số. Lần thử đầu tiên sử dụng các hệ số được hiển thị trong câu hỏi dưới dạng giá trị khởi đầu (ngoại trừ một giá trị) nhưng điều đó không thành công, do đó giảm C đã được thử và hệ số kết quả được đưa vào phù hợp thứ hai với C = 1.

st <- list(a = 43, b = -14, c = 25, B = 18) 
fm <- nls(y ~ (a + b * x + c * x^2)/(9 + B * x), start = st) 
fm2 <- nls(y ~ (a + b * x + c * x^2)/(9 + B * x + C * x^2), start = c(coef(fm), C = 1)) 

plot(y ~ x) 
lines(fitted(fm2) ~ x, col = "red") 

(tiếp tục sau khi biểu đồ)

Lưu ý: Dưới đây là một ví dụ của việc sử dụng nls2 để có được giá trị bắt đầu với việc tìm kiếm ngẫu nhiên. Chúng tôi giả định rằng các hệ số nằm giữa -50 và 50.

library(nls2) 

set.seed(123) # for reproducibility 
v <- c(a = 50, b = 50, c = 50, B = 50, C = 50) 
st0 <- as.data.frame(rbind(-v, v)) 
fm0 <- nls2(y ~ (a + b * x + c * x^2)/(9 + B * x + C * x^2), start = st0, 
    alg = "random", control = list(maxiter = 1000)) 

fm3 <- nls(y ~ (a + b * x + c * x^2)/(9 + B * x + C * x^2), st = coef(fm0)) 

Answer 2

Vì bạn đã có một dự đoán lý thuyết, bạn không dường như cần một mô hình mới và nó thực sự chỉ là một công việc vẽ sơ đồ:

png(); plot(y~x) 
lines(x,mod,col="blue") 
dev.off() 

Bạn không thể mong đợi lm để tạo ra một xấp xỉ tốt cho một vấn đề phi tuyến tính. Mẫu số liên quan đến x trong biểu thức lý thuyết đó làm cho phi tuyến vốn có.