-6

I have a time series data set and each time series has datapoint of 30-year from different/same species. I am developing a forecasting model using the first 23 years of data from each time series data point and I am using the rest 7 years as test set to know the predictive ability of model but the nonlinear model (Model 6 and Model 7) don't give succinct result?

Data:

DD <- structure(list(Plot = structure(c(1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 
3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 
3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 
3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 
4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 
4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L), .Label = c("A", 
"B", "C", "D"), class = "factor"), Species = structure(c(2L, 
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L), .Label = c("BD", "BG"), class = "factor"), Year = c(37L, 
38L, 39L, 40L, 41L, 42L, 43L, 44L, 45L, 46L, 47L, 48L, 49L, 50L, 
51L, 52L, 53L, 54L, 55L, 56L, 57L, 58L, 59L, 60L, 61L, 62L, 63L, 
64L, 65L, 66L, 37L, 38L, 39L, 40L, 41L, 42L, 43L, 44L, 45L, 46L, 
47L, 48L, 49L, 50L, 51L, 52L, 53L, 54L, 55L, 56L, 57L, 58L, 59L, 
60L, 61L, 62L, 63L, 64L, 65L, 66L, 37L, 38L, 39L, 40L, 41L, 42L, 
43L, 44L, 45L, 46L, 47L, 48L, 49L, 50L, 51L, 52L, 53L, 54L, 55L, 
56L, 57L, 58L, 59L, 60L, 61L, 62L, 63L, 64L, 65L, 66L, 37L, 38L, 
39L, 40L, 41L, 42L, 43L, 44L, 45L, 46L, 47L, 48L, 49L, 50L, 51L, 
52L, 53L, 54L, 55L, 56L, 57L, 58L, 59L, 60L, 61L, 62L, 63L, 64L, 
65L, 66L, 37L, 38L, 39L, 40L, 41L, 42L, 43L, 44L, 45L, 46L, 47L, 
48L, 49L, 50L, 51L, 52L, 53L, 54L, 55L, 56L, 57L, 58L, 59L, 60L, 
61L, 62L, 63L, 64L, 65L, 66L, 37L, 38L, 39L, 40L, 41L, 42L, 43L, 
44L, 45L, 46L, 47L, 48L, 49L, 50L, 51L, 52L, 53L, 54L, 55L, 56L, 
57L, 58L, 59L, 60L, 61L, 62L, 63L, 64L, 65L, 66L, 37L, 38L, 39L, 
40L, 41L, 42L, 43L, 44L, 45L, 46L, 47L, 48L, 49L, 50L, 51L, 52L, 
53L, 54L, 55L, 56L, 57L, 58L, 59L, 60L, 61L, 62L, 63L, 64L, 65L, 
66L), Count = c(81L, 45L, 96L, 44L, 24L, 8L, 28L, 32L, 39L, 29L, 
40L, 17L, 4L, 12L, 18L, 11L, 63L, 98L, 78L, 76L, 67L, 36L, 56L, 
43L, 81L, 8L, 14L, 20L, 25L, 19L, 135L, 91L, 171L, 88L, 59L, 
1L, 1L, 1L, 2L, 1L, 11L, 9L, 34L, 15L, 32L, 21L, 33L, 43L, 39L, 
20L, 6L, 3L, 9L, 9L, 28L, 16L, 15L, 2L, 1L, 1L, 34L, 16L, 19L, 
35L, 32L, 7L, 2L, 30L, 29L, 25L, 28L, 11L, 31L, 31L, 28L, 27L, 
34L, 110L, 87L, 103L, 72L, 19L, 46L, 43L, 107L, 32L, 26L, 31L, 
12L, 29L, 23L, 40L, 50L, 23L, 34L, 11L, 9L, 4L, 24L, 55L, 14L, 
16L, 51L, 43L, 2L, 13L, 8L, 96L, 52L, 118L, 32L, 1L, 8L, 17L, 
34L, 29L, 38L, 15L, 4L, 38L, 2L, 1L, 1L, 1L, 1L, 1L, 3L, 3L, 
4L, 6L, 4L, 4L, 10L, 6L, 7L, 9L, 15L, 30L, 25L, 36L, 13L, 17L, 
43L, 36L, 60L, 50L, 26L, 13L, 13L, 27L, 18L, 56L, 96L, 16L, 54L, 
2L, 2L, 9L, 5L, 5L, 6L, 2L, 6L, 2L, 3L, 4L, 3L, 136L, 71L, 116L, 
28L, 23L, 76L, 64L, 98L, 58L, 26L, 13L, 13L, 13L, 18L, 19L, 24L, 
18L, 17L, 3L, 23L, 19L, 9L, 11L, 13L, 20L, 29L, 29L, 17L, 20L, 
26L, 71L, 63L, 53L, 54L, 20L, 22L, 18L, 93L, 50L, 18L, 12L, 12L, 
31L), LogCount = c(1.908385019, 1.653212514, 1.982271233, 1.643462676, 
1.380211242, 0.903089987, 1.447158031, 1.505109978, 1.591064607, 
1.462397998, 1.602059991, 1.230448921, 0.602059991, 1.079181206, 
1.255272505, 1.041392685, 1.799340549, 1.991226076, 1.892094603, 
1.880813592, 1.826074803, 1.556302501, 1.748188027, 1.633468456, 
1.908485019, 0.903089987, 1.146128035, 1.301029996, 1.397940009, 
1.278753601, 2.130333768, 1.95904139, 2.2329961, 1.94448267, 
1.770852012, 0, 0, 0, 0.30102999, 0, 1.0411392685, 0.954242509, 
1.531478917, 1.176031259, 1.505149978, 1.322219295, 1.51851394, 
1.6334684456, 1.591064607, 1.301029996, 0.77815125, 0.477121255, 
0.954242509, 0.954242509, 1.447158031, 1.204119983, 1.176091259, 
0.301029996, 0, 0, 1.531478917, 1.204119983, 1.278753501, 1.544068044, 
1.505149978, 0.084509804, 0.301029996, 1.477121255, 1.462397998, 
1.397940009, 1.447158031, 1.041392685, 1.491361694, 1.491361694, 
1.447158031, 1.431363754, 1.531478917, 2.041392685, 1.939519253, 
2.012837225, 1.857332495, 1.278753601, 1.662757382, 1.633468456, 
2.029383778, 1.505149978, 1.414973348, 1.491361594, 1.079181245, 
1.462397998, 1.361727835, 1.602059991, 1.698970004, 1.361727836, 
1.531478917, 1.041392685, 0.954242509, 0.602059991, 1.380211242, 
1.740362689, 1.146128036, 1.204119983, 1.707570176, 1.633468456, 
0.301029996, 1.113943352, 0.903089987, 1.982271233, 1.716003344, 
2.071882007, 1.50514997, 0, 0.903089987, 1.230448921, 1.53147891, 
1.2397998, 1.57978359, 1.176091259, 0.602059991, 1.57978359, 
0.301029996, 0, 0, 0, 0, 0, 0.477121255, 0.477121255, 0.602059991, 
0.77815125, 0.602059991, 0.602059991, 1, 0.77815125, 0.84509804, 
0.95424509, 1.176091259, 1.477121255, 1.39790009, 1.555302501, 
1.113943352, 1.230448921, 1.633468456, 1.555302501, 1.77815125, 
1.698970004, 1.414973348, 1.113943352, 1.113943352, 1.431353754, 
1.255272505, 1.748188027, 1.982271233, 1.204119983, 1.73239376, 
1.431363754, 1.361727835, 0.954242509, 0.698970004, 0.698970004, 
0.77815125, 0.301029996, 0.77815125, 0.301029996, 0.477121255, 
0.602059991, 0.477121255, 2.133538908, 1.851258349, 2.064457989, 
1.447158031, 1.361727836, 1.880813592, 1.806179974, 1.991226076, 
1.763427994, 1.414973348, 1.113943352, 1.113943352, 1.113943352, 
1.255272505, 1.278753601, 1.380211242, 1.255272505, 1.230446921, 
0.477121255, 1.361727835, 1.278753601, 0.954242509, 1.0411392685, 
1.113943352, 1.301029996, 1.462397998, 1.462397998, 1.230448921, 
1.301029995, 1.414973348, 1.851258349, 1.799340549, 1.72427587, 
1.73239376, 1.301029996, 1.342422681, 1.255272505, 1.968482949, 
1.698970004, 1.255272505, 1.079181246, 1.079181246, 1.491361694
), Diff = c(-0.255272505, 0.329058719, -0.338818557, -0.263241434, 
-0.077121255, 0.544068044, 0.057991947, 0.085910629, -0.128666609, 
0.139661993, -0.37161107, -0.62838893, 0.477121255, 0.176091259, 
-0.21387982, 0.757947864, 0.191885527, -0.099131473, -0.011281011, 
-0.054738789, -0.269772302, 0.191885526, -0.114719571, 0.275016563, 
-1.005395032, 0.243038049, 0.15490196, 0.096910013, -0.119186408, 
NA, -0.171292376, 0.273954718, -0.288513438, -0.17363066, -1.770852012, 
0, 0, 0.301029996, -0.301029996, 1.041392685, -0.087150176, 0.577235408, 
-0.355387658, 0.329058719, -0.182930683, 0.196294545, 0.110954516, 
-0.042403849, -0.290034611, -0.522878746, -0.301029995, 0.477121254, 
0, 0.492915522, -0.243038048, -0.028028724, -0.875061263, -0.301029996, 
0, 1.531078917, -0.32735893, 0.070633618, 0.265310043, -0.038918066, 
-0.660051938, -0.544068044, 1.176091259, -0.014723257, -0.064457989, 
0.049218022, -0.405765346, 0.449969009, 0, -0.044203663, -0.015794267, 
0.100115153, 0.509913768, -0.101873432, 0.073317972, -0.155504729, 
-0.578578895, 0.384054231, -0.029289376, 0.395915322, -0.5202338, 
-0.09017663, 0.076388346, -0.412180448, 0.383216752, -0.100670162, 
0.240332155, 0.096910013, -0.337242168, 0.169751081, -0.490086232, 
-0.087150176, -0.352182518, 0.778151251, 0.360151447, -0.594234653, 
0.057991947, 0.503450193, -0.07410172, -1.33243846, 0.812913356, 
-0.210853365, 1.079181246, -0.266267889, 0.355878663, -0.566732029, 
-1.505149978, 0.903089987, 0.327358934, 0.301029996, -0.069080919, 
0.117385599, -0.403692338, -0.574031268, 0.977723606, -1.278753601, 
-0.301029996, 0, 0, 0, 0, 0.477121255, 0, 0.124938736, 0.176091259, 
-0.176091259, 0, 0.397490009, -0.2218485, 0.06690679, 0.10914469, 
0.22184875, 0.301029996, -0.079181206, 0.158362092, -0.442359149, 
0.116505569, 0.403019535, -0.077165955, 0.221848749, -0.079181206, 
-0.283996656, -0.301029996, 0, 0.317420412, -0.176091259, 0.492915522, 
0.23483206, -0.77815125, 0.528273777, -0.301029996, -0.069635928, 
-0.407485327, -0.255272505, 0, 0.079181246, -0.477121254, 0.477121254, 
-0.477121254, 0.176091259, 0.124938736, -0.124938736, 1.656417653, 
-0.282280559, 0.21319964, -0.617299958, -0.085430195, 0.5191085756, 
-0.074533518, 0.185045102, -0.227798082, -0.348454546, -0.301029996, 
0, 0, 0.141329153, 0.023481096, 0.101457641, -0.124938737, -0.024823584, 
-0.753327666, 0.884606581, -0.082974235, -0.324511092, 0.087150176, 
0.072550667, 0.187086644, 0.161368002, 0, -0.231949077, 0.070581075, 
0.113903352, 0.436285001, -0.00519178, -0.075054679, 0.00811789, 
-0.431363764, 0.041392685, -0.087150176, 0.713210444, -0.269512945, 
-0.443697499, -0.176091259, 0, 0.412180448, -0.148939013)), class = "data.frame", row.names = c(NA, 
-210L))

Code:

for(f in 1:11){
  for(b in 1:5){
    for (c in 1:5){

      #To select test sets 1,2,3,4, and 5 years beyond the training set:

      #Calculate the mean of abundance for the training set years.
      Model1<-lm(mean~1, data=DD1)

      #

Output2:

2 3 0.676209994477288 1.9365051784348e-09 4.44089209850063e-16
3 53 11.9236453578109 2.06371097988267e-09 1.13686837721616e-13
4 31 1.94583877614293 1.11022302462516e-15 1.99840144432528e-15
5 4 8.06660449042397 1.48071350736245e-08 3.19744231092045e-14
6 5 10.5321102149558 9.31706267692789e-10 1.4210854715202e-14

..
Stackuser
  • 59
  • 7
  • 5
    There are a few important steps you can take to improve the quality of this post (which will, in turn, get you better answers, faster). (1) include example data that can be copied, inline, as code, not as a screenshot. (2) include expected output - that helps people check their answers against your expectation. (3) ask a more specific question than "better code, please?" - try to narrow down your question to a specific issue you're having trouble with. See [mcve] for more. – andrew_reece May 03 '19 at 22:00
  • Could you execute `dput(DD)` and add the output of this function into your question? – Artem May 04 '19 at 07:46
  • It looks somewhere in your algorithm `nls` function is fed with one-row data.frame. So `nls` throws "matrix is singular" . It seems on some steps you are trying to fit the data with only one point. – Artem May 06 '19 at 05:42
  • Can you please help me in fixing what you have seen, Artem? – Stackuser May 06 '19 at 11:17

1 Answers1

0

First, please see the time series graph of counts for different species and plots below.

library(ggplot2)
ggplot(DD, aes(Year, Count)) +
  geom_point() +
  geom_line() +
  facet_grid(Plot ~ Species) + 
  scale_y_log10()

Count vs year

It is seen that there is no obvious trend which can be fitted by power or log-power function using nls.

Second, as I understand you are trying to use nls to predict outside the training data set. Usually it is not quite an effective to use least square models because of auto-correlated nature of time-series.

Third, the simplest prediction algorithm is Holt-Winters (see "dirty" implementation below). You can use as well a ton of other algorithms like ARIMA, exponential smoothing state space etc.

x <- ts(DD[DD$Species == "BG" & DD$Plot == "elq1a3", ]$LogCount)
m <- HoltWinters(x, gamma = FALSE)
library(forecast)
f <- forecast(m, 2)
plot(f, main = "elq1a3 at BG")

HoltWinters

Fourth, about your algorithm in question, it throws

Error in qr.solve(QR.B, cc) : singular matrix 'a' in solve.

The reason is that in the first step of for-loop (at f = b = c = 1 DD2 data frame contains just one row. And executing

Model6<-nls(Diff~1+Count^T,start=list(T=1),trace=TRUE,algorithm ="plinear",data=DD2)

means that you are trying to fit a curve using only one data point, which is impossible. However if you change f value in for-loop from 1:11 to 2:11 another error thrown:

Error in nls(Diff ~ 1 + Count^T, start = list(T = 1), trace = TRUE, algorithm = "plinear", : step factor 0.000488281 reduced below minFactor 0.000976562.

In this case you cannot use "naive" approach used by plinear algorithm with self-starting inital value and, e.g. nls.control(min.factor = 1e-5). You must feed all initial coefficients explicitely with default Gauss-Newton algorithm. Quite exausting, please try yourself :)

Artem
  • 3,304
  • 3
  • 18
  • 41