CODEWITHSUNDEEP
rlogistic-regressionpanel-datanls
Lovely97
Lovely97

Reputation: 1

R Finding logistic curve with nls

I have a problem with logistic curve in R for panel data which are here:

https://docs.google.com/spreadsheets/d/1SO3EzFib7T3XqTz1xZCU2bZFTB0Ddhou/edit?usp=sharing&ouid=110784858039906954607&rtpof=true&sd=true

I tried:

log <- nls(lrc~SSlogis(time, Asym, xmid, scal), data = data_log)

I have an error: 'qr.solve(QR.B, cc)':singular matrix 'a' in solve. What can I do?

Upvotes: 0

Views: 293

Answers (1)

Ben Bolker
Ben Bolker

Reputation: 226692

I get a different error:

Error in nls(y ~ 1/(1 + exp((xmid - x)/scal)), data = xy, start = list(xmid = aux[[1L]], : step factor 0.000488281 reduced below 'minFactor' of 0.000976562

(it's not surprising that different platforms will get slightly different results on numerically "difficult" problems ...)

However, here's what plot(lrc ~ time, data = dd) produces when I use your data:

plot of lrc vs time - no pattern

It seems optimistic to think that you could fit a logistic curve to these data, or that the fit would make very much sense ...

I did find that I could fit a logistic to the logged data, i.e. nls(log(lrc) ...)

plot(log(lrc) ~ time, data =dd)
tvec <- seq(2008, 2016, by = 0.1)
lines(tvec, predict(m2, newdata=data.frame(time=tvec)), col=2, lwd=2)

log scale plot + prediction

If I really needed logistic coefficients for this plot (e.g. to compare with other cases) I would fit a linear model to the data, assume that the midpoint was equal to the midpoint of the data, set the scaling parameter equal to the linear slope coefficient divided by 4 (this is a standard rule of thumb for the logistic), and say that the asymptotic value could not be estimated.

Plotting all of your data unit-by-unit doesn't make much more hopeful that you will be able to fit logistic curves unit-by-unit either. While there might be a few units where a logistic curve is a sensible description of the data, it's not in most cases. You might have to back up and consider your analytical strategy — i.e., what are you hoping to learn from these data, and how might you go about it? (If you can frame the question suitably you could post on CrossValidated. If you are a student/trainee, you might want to ask your supervisor/teacher/mentor for advice ...

panel plot: all units

library(readxl)
library(tidyverse)
(dd <- read_excel("data2.xlsx")
  %>% pivot_longer(-code, names_transform = as.numeric,
                   names_to = "year")
  ## indices to break groups into chunks/facets
  %>% mutate(grp_cat = factor(as.numeric(factor(code)) %% 30))
)
gg1 <- ggplot(dd, aes(year, value, group = code,
               colour=code)) + geom_point() +
  facet_wrap(~grp_cat, scale="free_y") +
  expand_limits(y=0) +
  theme_bw() +
  theme(legend.position = "none",
        panel.spacing = grid::unit(0, "lines"),
        axis.text.x = element_blank())
gg1 + geom_smooth(se=FALSE)

ggsave("all.png", width=8, height=8)

Upvotes: 1

Related Questions