Manuel Popp
Manuel Popp

Reputation: 1175

Correct use of spaMM::Matern in spatial GLMM

I want to fit a GLM and control for spatial relationships within the data. To this end, I wanted to include the term spaMM::Matern(1 | x + y) where x and y are the column and row index from the raster stack where the other data points were sampled.

Without the term, I can sucessfully fit a model using spaMM::fitme(); however, once I introduce the Matérn covariance term

mod_matern <- spaMM::fitme(
  response ~ type * (
    cmi_dry_rate + cmi_mean + pr_mean + swb_mean + tas_max + tas_mean + tas_min
  ) + spaMM::Matern(1 | x + y),
  data = df, family = spaMM::negbin()
)

running the code fails with the error:

Error in if (term1 == "multIMRF") { : the condition has length > 1

I tried checking for duplicates in the coordinates or NA values, etc. with no success.

Am I using the function incorrectly, is there a bug, or some other issue?

Here is a MWE:

library(spaMM)

df <- structure(
  list(
    response = c(0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0),
    cmi_dry_rate = c(0.303921580314636, 0.220588237047195, 0.426470577716827,
     0.220588237047195, 0.0833333358168602, 
     0.549019634723663, 0.338235288858414, 0.343137264251709, 0.485294103622437, 
     0.519607841968536, 0.460784316062927, 0.480392158031464, 0.279411762952805, 
     0.333333343267441, 0.480392158031464, 0.470588237047195, 0.333333343267441, 
     0.436274498701096, 0.397058814764023, 0.147058829665184),
    cmi_mean = c(15.6495094299316, 
                35.3279418945312, -1.75882351398468, 33.1259803771973, 82.2848052978516, 
                -15.6254901885986, 16.6882362365723, 8.32058811187744, -5.65882349014282, 
                -11.6999998092651, -5.71029424667358, -8.84951019287109, 22.1514701843262, 
                7.43235301971436, -6.78333330154419, -5.53725481033325, 5.11617660522461, 
                -4.68088245391846, -0.450490206480026, 65.5897064208984),
    pr_mean = c(5596.2890625, 
                8757.5732421875, 4597.5537109375, 8431.671875, 14661.2353515625, 
                4254.0390625, 6983.82373046875, 4515.89208984375, 4661.51953125, 
                4162.5390625, 4436.529296875, 3783.81372070312, 6743.2255859375, 
                5202.90185546875, 4468.42138671875, 4487.48046875, 4415.81884765625, 
                4606.16650390625, 5045.90673828125, 10231.9267578125),
    swb_mean = c(124, 
                164.470581054688, -35.235294342041, 176.882354736328, 193.058822631836, 
                -190.529418945312, 129.235290527344, 71.529411315918, -83, -138.823532104492, 
                -79.8823547363281, -99.529411315918, 132.117645263672, 41.8823547363281, 
                -85.3529434204102, -70.8823547363281, 35.764705657959, -59.294116973877, 
                -73.2352905273438, 189.882354736328),
    tas_max = c(2934, 2908, 
                2943, 2922, 2911, 2907, 2924, 2923, 2956, 2948, 2943, 2917, 2934, 
                2934, 2940, 2931, 2908, 2938, 2903, 2919),
    tas_mean = c(2745.087890625, 
                2734.55102539062, 2734.37963867188, 2745.912109375, 2785.23608398438, 
                2771.5185546875, 2740.6806640625, 2644.875, 2711.95825195312, 
                2738.53247070312, 2737.27783203125, 2720.97680664062, 2757.03247070312, 
                2698.49072265625, 2729.6943359375, 2709.12036132812, 2666.7685546875, 
                2736.75, 2703.27783203125, 2707.5),
    tas_min = c(2508, 2500, 2401, 
                2517, 2664, 2568, 2494, 2287, 2431, 2419, 2418, 2487, 2560, 2355, 
                2467, 2436, 2418, 2461, 2430, 2458),
    type = structure(
      c(2L, 1L, 2L, 1L, 2L, 1L, 1L, 2L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 2L, 2L),
      levels = c("loss", "stable"), class = "factor"),
    y = c(5901, 8724, 7145, 9156, 9423, 8102, 8765, 5602, 8162, 7469, 7350, 6686, 
          5486, 6067, 7960, 7958, 5169, 7636, 7912, 5922),
    x = c(50732, 24826, 57844, 24456, 27604, 12517, 20260, 62334, 68779, 57841, 
         57326, 13193, 45963, 57225, 62833, 63908, 9151, 62284, 20826, 
         53078)
    ),
  row.names = c(9595L, 158L, 7247L, 3822L, 5852L, 4028L, 4672L, 5131L, 4826L,
                4352L, 3341L, 5385L, 5619L, 9913L, 5182L, 
                8859L, 9261L, 4770L, 9140L, 7469L), class = "data.frame"
  )

mod_matern <- spaMM::fitme(
  response ~ type * (
    cmi_dry_rate + cmi_mean + pr_mean + swb_mean + tas_max + tas_mean + tas_min
  ) + spaMM::MaternCorr(1 | x + y),
  data = df, family = spaMM::negbin()
)

Upvotes: 0

Views: 35

Answers (1)

Manuel Popp
Manuel Popp

Reputation: 1175

As @user206550 suggested, using just Matern(1 | x + y) works. It seems strange that spaMM::Matern(1 | x + y) would cause the kind of error message mentioned, but it seems it just is like that.

Upvotes: 2

Related Questions