Ошибка в коде для логистической регрессии

Question

Я пишу следующий сценарий, используя PopulationGrowth.csv:

dat <-read.csv("/Path/PopulationGrowth.csv") 
class = "data.frame", row.names=c(NA,-16L) 
m1 <- nls(pop~SSlogis(Year,asym,xmid,scal),data=dat) 
par(las=1,bty="l",mar=c(5,6,2,2)+0.1)  ## graphics tweaks 
with(dat,plot(CentralOakland~Year,ylab="")) 
mtext("Population",side=2,las=0,line=4) 
yearvec <- 1940:2010 
lines(yearvec,predict(m1,newdata=data.frame(Year=yearvec))) 

После последней строки я получаю следующую ошибку:

Error in predict(m1, newdata = data.frame(Year = yearvec)) : object 'm1' not found 

Вот данные я работа с:

Year CentralOakland 
1940 7852 
1950 8452 
1960 6701 
1970 6135 
1980 5872 
1990 5406 
2000 5281 
2010 6086 

Я также пытаюсь прогнозировать тенденции народонаселения в течение следующих 30 лет. Есть основания полагать, что население будет расти в течение следующих 30 лет. Какую функцию необходимо предсказать?

медведь со мной, я новичок в R.

Answer 1

Есть несколько маленьких глюков здесь:

• вы не дали нам воспроизводимый пример: смотри, например, http://tinyurl.com/reproducible-000 ...
• вы, вероятно, хотите read.csv, а не read.table
• Я не знаю, что quartz=(title...) делает. Может быть quartz(title=...)? В любом случае, с помощью plot() автоматически открыть новое графическое окно
• использованием attach часто плохая идея
• использование lines вместо curve, чтобы добавить новый (например) прогнозируемые данные, установленные на участке. Для конкретного случая линейной регрессии (см. Следующий пункт) вы можете использовать abline.
• на сегодняшний день самое главное: то, что вы делаете здесь не логистической регрессии - он использует glm, чтобы соответствовать обычной линейной модели (регрессии по методу наименьших квадратов), который также может быть достигнуто более легко с помощью lm(). Логистическая регрессия использовала бы family=binomial, а не family=gaussian. Кроме того, логистическая регрессия предназначена для переменных бинарного ответа, так это вряд ли целесообразна для моделирования населения Центральной Окленде ... Я подозреваю, что вы путаете логистическую регрессию и логистический кривой рост ...

Вот некоторые пересмотренный код:

dat <- read.csv("/Path/PopulationGrowth.csv") 
with(dat,plot(Year,CentralOakland,xlab="Year",ylab="Population", 
    main="Central Oakland Population vs. Year")) 
g <- lm(CentralOakland~Year,data=dat) 
abline(g) 
## OR 
yearvec <- seq(min(Year),max(Year),length=51) 
lines(yearvec,predict(g,newdata=data.frame(Year=yearvec),type="response")) 

с другой стороны, может быть, это то, что вы действительно хотели сделать (вы не дали нам данные, поэтому я пошел в Википедию и получил некоторые ...)

## from http://en.wikipedia.org/wiki/Oakland,_California 
dat <- structure(list(year = c(1860L, 1870L, 1880L, 1890L, 1900L, 1910L, 
1920L, 1930L, 1940L, 1950L, 1960L, 1970L, 1980L, 1990L, 2000L, 
2010L), pop = c(1543L, 10500L, 34555L, 48682L, 66960L, 150174L, 
216261L, 284063L, 302163L, 384575L, 367548L, 361561L, 339337L, 
372242L, 399484L, 390724L)), .Names = c("year", "pop"), 
    class = "data.frame", row.names = c(NA,-16L)) 

Место:

## SSlogis() is a "self-starting" logistic 
m1 <- nls(pop~SSlogis(year,asym,xmid,scal), 
      data=dat) 

Черчение:

par(las=1,bty="l",mar=c(5,6,2,2)+0.1) ## graphics tweaks 
with(dat,plot(pop~year,ylab="")) 
## add y-label separately 
mtext("population",side=2,las=0,line=4) 
yearvec <- 1860:2010 
lines(yearvec,predict(m1,newdata=data.frame(year=yearvec)))