Below is what I have done and the error message I am getting each time

gee_fit <- geeglm(rslt~ ns(day,3), data = epidem,id =patid, family = binomial, corstr = 'exchangeable')

To obtain predicted probability for selected day values

predict(gee_fit, newdata = data.frame(day = c(20, 30)),  type="response")

If I use the below to get the standard deviation, I get an error message

day_grid<-seq(from=min(day), to=max(day))
preds<-predict(gee_fit, newdata=list(day_grid), se=T)

Error in XRinv^2 %*% rep(res.var, p) :
requires numeric/complex matrix/vector arguments

Where am I going wrong?

Data

structure(list(day = structure(c(2, 8, 13, 15, 12, 29, 33, 22, 
17, 21, 44, 18, 22, 26, 8, 13, 19, 5, 10, 16, 41, 21, 26, 21, 
2, 20, 28, 34, 29, 30, 37, 14, 20, 21, 12, 19, 8, 9, 14, 19, 
14, 18, 12, 17, 24, 12, 16, 12, 33, 60, 42, 48, 52, 58, 51, 57, 
64, 24, 37, 44, 19, 25, 5, 11, 33, 17, 22, 28, 27, 32, 38, 2, 
7, 12, 19, 35, 41, 50, 3, 14, 59, 47, 53, 61, 37, 17, 8, 13, 
21, 13, 5, 13, 2, 8, 14, 56, 84, 10, 32, 11, 7, 11, 10, 57, 11, 
11, 12, 27, 9, 13, 10, 15, 11, 15, 9, 12, 17, 7, 11, 13, 7, 27, 
38, 15, 17, 10, 12, 19, 4, 13, 18, 11, 65, 15, 26, 2, 22, 28, 
28, 40, 60, 74, 9, 12, 8, 13, 18, 7, 65, 9, 13, 18, 14, 8, 12, 
17, 20, 25, 19, 8, 10, 12, 7, 35, 56, 5, 7, 11, 14, 5, 16, 7, 
12, 19, 24, 41, 2, 26, 2, 4, 13, 15, 12, 17, 8, 6, 11, 17, 9, 
21, 30, 35, 8, 15, 4, 25, 13, 47, 38, 23, 14, 16, 21, 7, 12, 
18, 9, 13, 24, 25, 33, 41, 7, 16, 34, 50, 53, 49, 47, 41, 43, 
54, 68, 32, 50, 36, 54, 11, 17, 11, 16, 35, 10, 12, 26, 43, 6, 
12, 15, 6, 10, 37, 46, 57, 59, 53, 42, 55, 53, 7, 13, 20, 7, 
10, 13, 57, 50, 56, 41, 43, 52, 41, 55, 40, 56, 41, 46, 74, 21, 
35, 46, 58, 52, 61, 58, 12, 17, 20, 25, 51, 56, 63, 24, 29, 35, 
25, 46, 60, 45, 61, 63, 72, 65, 73, 63, 61, 53, 78, 11, 30, 23, 
8, 13, 2, 9, 44, 75, 38, 41, 54, 57, 24, 60, 20, 26, 7), label = "day", format.stata = "%9.0g"), 
    rslt = structure(c(0, 0, 1, 1, 0, 0, 0, 1, 0, 1, 1, 
    0, 1, 1, 0, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 
    1, 0, 0, 1, 1, 1, 0, 0, 1, 1, 1, 1, 0, 1, 1, 1, 1, 0, 1, 
    1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
    1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 0, 1, 1, 0, 1, 1, 1, 1, 0, 
    0, 1, 0, 0, 1, 0, 1, 0, 1, 0, 1, 1, 1, 0, 0, 1, 1, 1, 1, 
    1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 
    0, 1, 1, 0, 0, 0, 1, 1, 1, 1, 0, 0, 0, 1, 1, 1, 1, 1, 0, 
    1, 1, 1, 0, 0, 1, 1, 1, 0, 0, 0, 1, 1, 1, 1, 0, 0, 1, 0, 
    1, 1, 1, 1, 1, 1, 0, 1, 0, 1, 1, 1, 0, 0, 0, 0, 0, 0, 1, 
    1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 0, 1, 0, 1, 1, 
    1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 
    1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 
    0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 1, 1, 1, 1, 
    1, 1, 1, 1, 1, 1, 1, 0, 0, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 
    1, 1, 1, 1, 1, 0, 0, 0, 1, 1, 1, 1, 1, 1, 0, 0, 1, 0, 1, 
    1, 0, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 
    1), label = "test rslts", format.stata = "%9.0g"), 
    patid = structure(c(147, 147, 147, 148, 149, 150, 150, 151, 
    152, 152, 152, 153, 153, 153, 154, 154, 154, 155, 155, 155, 
    156, 157, 157, 158, 159, 160, 160, 160, 161, 162, 162, 163, 
    163, 164, 165, 166, 167, 168, 169, 169, 170, 170, 171, 171, 
    171, 172, 172, 173, 174, 175, 176, 176, 177, 177, 178, 178, 
    178, 179, 179, 179, 181, 181, 182, 182, 183, 184, 184, 184, 
    185, 185, 185, 186, 186, 186, 187, 188, 188, 188, 189, 189, 
    190, 191, 191, 191, 192, 193, 194, 194, 194, 195, 196, 197, 
    198, 198, 198, 199, 200, 72, 72, 16, 77, 83, 17, 17, 108, 
    84, 18, 18, 81, 81, 71, 66, 111, 111, 140, 20, 4, 21, 22, 
    22, 107, 139, 55, 19, 116, 23, 23, 23, 13, 13, 13, 117, 117, 
    82, 82, 104, 104, 104, 96, 96, 96, 96, 138, 110, 86, 86, 
    86, 5, 5, 114, 114, 114, 130, 3, 3, 3, 3, 3, 67, 106, 106, 
    106, 65, 65, 24, 12, 12, 12, 12, 85, 133, 15, 15, 15, 91, 
    1, 144, 144, 57, 57, 57, 105, 100, 100, 56, 11, 11, 11, 124, 
    60, 60, 60, 14, 14, 10, 52, 9, 9, 103, 68, 128, 29, 29, 54, 
    54, 54, 61, 61, 143, 27, 27, 27, 28, 30, 95, 123, 33, 73, 
    137, 134, 97, 97, 97, 118, 119, 120, 35, 31, 31, 34, 34, 
    34, 32, 126, 38, 38, 59, 59, 59, 76, 76, 94, 78, 36, 109, 
    37, 93, 115, 141, 99, 99, 99, 75, 75, 75, 62, 39, 113, 64, 
    146, 121, 92, 80, 40, 87, 25, 25, 25, 102, 102, 44, 45, 46, 
    74, 136, 131, 131, 63, 63, 89, 89, 89, 42, 42, 42, 127, 58, 
    48, 47, 43, 145, 125, 2, 122, 129, 8, 90, 135, 112, 49, 41, 
    142, 142, 53, 53, 69, 98, 101, 132, 88, 50, 7, 26, 51, 6, 
    70), label = "group(record)", format.stata = "%9.0g")), row.names = c(NA, 
-316L), class = c("tbl_df", "tbl", "data.frame"))