Как найти кривую, которая соответствует серии точек на R? - PullRequest
Новая
       94

Как найти кривую, которая соответствует серии точек на R?

0 голосов
/ 26 марта 2020

Мне нужно выяснить уравнение кривой мощности, которая приспосабливается к загрязнению в день определенной болезни, чтобы я мог сделать прогноз, данные следующие: 

Day     Contaminated

26/feb  1
29/feb  2
04/mar  3
05/mar  8
06/mar  13
07/mar  19
08/mar  25
10/mar  34
11/mar  52
12/mar  81
13/mar  98
14/mar  121
15/mar  176
16/mar  234
17/mar  291
18/mar  428
19/mar  621
20/mar  904
21/mar  1128
22/mar  1546
23/mar  1891
24/mar  2201
25/mar  2433

Я считаю, что я нужно сделать регрессию кривой мощности (NonLinearRegression) в R, но я не знаю, как это реализовать.

1 Ответ

0 голосов
/ 27 марта 2020

Вот подход, использующий data.table, ggplot2 и nls.

Сначала давайте исправим даты в стандартном формате и преобразуем их в целые числа, чтобы мы могли выполнить некоторые вычисления. 

library(data.table)
library(ggplot2)
setDT(data)
data[,Day:= as.Date(Day,"%d/%b")]
data[,Int := as.integer(Day)-min(as.integer(Day))]

Затем мы используем nls для подгонки модели к данным. Мы будем использовать формулу y = a * x ^ b. 

nls(formula = Contaminated ~ a * Int ^ b, data,start=list(a=1,b=1))
# Nonlinear regression model
#  model: Contaminated ~ a * Int^b
#   data: data
#        a         b 
#2.272e-05 5.571e+00 
# residual sum-of-squares: 123279
#
#Number of iterations to convergence: 48 
#Achieved convergence tolerance: 7.832e-07

Теперь мы можем просматривать результаты с помощью ggplot. 

ggplot(data, aes(x=Int,y=Contaminated)) + 
  geom_point() +
  scale_x_continuous(breaks = c(0,10,20), labels = data$Day[data$Int %in% c(0,10,20)]) + xlab("Date") +
  geom_smooth(method="nls", formula = y ~ a * x ^ b,method.args = list(start = c(a=1, b=1)),se=FALSE, linetype = 1)

enter image description here Данные 

data <- structure(list(Day = c("26/feb", "29/feb", "04/mar", "05/mar", 
"06/mar", "07/mar", "08/mar", "10/mar", "11/mar", "12/mar", "13/mar", 
"14/mar", "15/mar", "16/mar", "17/mar", "18/mar", "19/mar", "20/mar", 
"21/mar", "22/mar", "23/mar", "24/mar", "25/mar"), Contaminated = c(1L, 
2L, 3L, 8L, 13L, 19L, 25L, 34L, 52L, 81L, 98L, 121L, 176L, 234L, 
291L, 428L, 621L, 904L, 1128L, 1546L, 1891L, 2201L, 2433L)), class = "data.frame", row.names = c(NA, 
-23L))
...