Need GAM coefficients orderly within graphs in a grid

Question

I am running GAMs for a dataset that contains trends in 16 populations, and show the data and GAM trend lines for all 16 in a 4x4 grid. I need to include three of the GAM coefficients within each graph and haven't been able to come up with the correct script.

My df:

> dput(bb)
structure(list(Year = c(1970, 1971, 1972, 1973, 1974, 1975, 1976, 
1977, 1978, 1979, 1980, 1981, 1982, 1983, 1984, 1985, 1986, 1987, 
1988, 1989, 1990, 1991, 1992, 1993, 1994, 1995, 1996, 1997, 1998, 
1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009, 
2010, 2011, 2012, 2013, 2014, 2015, 2016, 2017, 2018), F = c(28, 
31, 27, 34, 24, 16, 18, 33, 31, 25, 33, 34, 46, 45, 45, 40, 18, 
18, 25, 14, 21, 16, 21, 20, 21, 20, 18, 24, 26, 17, 24, 16, 18, 
15, 21, 35, 53, 64, 51, 52, 59, 52, 56, 52, 50, 51, 44, 20, 30
), M = c(24, 39, 53, 39, 41, 18, 17, 33, 45, 48, 48, 60, 47, 
34, 36, 35, 28, 27, 24, 32, 29, 30, 35, 35, 38, 40, 37, 34, 43, 
31, 34, 44, 45, 46, 51, 61, 112, 109, 116, 119, 89, 106, 103, 
82, 87, 81, 67, 44, 35), U = c(7, 4, 3, 8, 7, 42, 29, 6, 21, 
14, 18, 7, 13, 13, 12, 40, 12, 51, 24, 0, 0, 0, 14, 79, 60, 64, 
67, 81, 96, 98, 68, 97, 118, 206, 156, 136, 1, 0, 0, 1, 0, 1, 
5, 1, 1, 4, 1, 71, 82), Tot = c(59, 74, 83, 81, 72, 76, 64, 72, 
97, 87, 99, 101, 106, 92, 93, 115, 58, 96, 73, 46, 50, 46, 70, 
134, 119, 124, 122, 139, 165, 146, 126, 157, 181, 267, 228, 232, 
166, 173, 167, 172, 148, 159, 164, 135, 138, 136, 112, 135, 147
), ratio = c(0.857142857142857, 1.25806451612903, 1.96296296296296, 
1.14705882352941, 1.70833333333333, 1.125, 0.944444444444444, 
1, 1.45161290322581, 1.92, 1.45454545454545, 1.76470588235294, 
1.02173913043478, 0.755555555555556, 0.8, 0.875, 1.55555555555556, 
1.5, 0.96, 2.28571428571429, 1.38095238095238, 1.875, 1.66666666666667, 
1.75, 1.80952380952381, 2, 2.05555555555556, 1.41666666666667, 
1.65384615384615, 1.82352941176471, 1.41666666666667, 2.75, 2.5, 
3.06666666666667, 2.42857142857143, 1.74285714285714, 2.11320754716981, 
1.703125, 2.27450980392157, 2.28846153846154, 1.50847457627119, 
2.03846153846154, 1.83928571428571, 1.57692307692308, 1.74, 1.58823529411765, 
1.52272727272727, 2.2, 1.16666666666667), popsize = c(NA, NA, 
NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
NA, NA, NA, NA, NA, NA, NA, NA, NA, 2074, 2074, 2074, 2074, 2074, 
2074, 2074, 1546, 1546, 1546, 1546, 1546, 1546, 1546, 1546, 1546, 
2826, 2826, 2826, 2826, 2826, 2826), rate = c(NA, NA, NA, NA, 
NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
NA, NA, NA, NA, NA, NA, NA, 6.70202507232401, 7.9556412729026, 
7.03953712632594, 6.07521697203472, 7.56991321118611, 8.72709739633558, 
12.8736740597878, 14.7477360931436, 15.006468305304, 10.7373868046572, 
11.1901681759379, 10.8020698576973, 11.1254851228978, 9.57309184993532, 
10.2846054333765, 10.608020698577, 4.77707006369427, 4.88322717622081, 
4.81245576786978, 3.96319886765747, 4.77707006369427, 5.20169851380042
)), row.names = c(NA, -49L), class = c("tbl_df", "tbl", "data.frame"
))

The script:

library(mgcv)

bb$rate <- (bb$Tot/bb$popsize)*100

m.bb <- gam(Tot~s(Year,k=6),family=poisson,data=bb)
mbb <- summary(m.bb)

YearP=seq(1970,2018,by=1)
mbb.pred=predict(m.bb,newdata=data.frame(Year=YearP),type="response",se.fit=T)

par(cex.lab=1.3,cex.axis=1.3)

plot(YearP, mbb.pred$fit, type="l", ylim=c(0,300),xlim=c(1970,2020),
     xlab="", ylab=" ",lwd=1.6,las=1)
points(bb$Year,bb$Tot,col=ifelse(bb$rate<4.5|is.na(bb$rate),'black','red'),pch=19,cex=1.7)
lines(YearP,mbb.pred$fit+2*mbb.pred$se.fit,lty="dotted",col="black")
lines(YearP,mbb.pred$fit-2*mbb.pred$se.fit,lty="dotted",col="black")
mtext('BB', side=3, line=0.8, at=1969,adj=0,cex=1)

r2 <- mbb$r.sq
rp = vector('expression',2)
rp[1] = substitute(expression(italic(R)^2 == r2), 
                   list(r2 = format(r2,dig=2)))[2]
legend('topright',adj=c(0.2,0.2),legend=rp,bty='n',cex=1.5)

The first part of the script runs the model, the second part plots the data, and the last part prints the edf for the GAM within the graph. But, I need both the edf, the adjusted R squared and the P-value within the graph, orderly on three lines above each other. I have copy-pasted the last part of the script from elsewhere, and I am unable to edit it the way I need.

Any help appreciated.

user_190277 · Answer

This is fussier than I would like but seems OK:

## extract values
r2 <- mbb$r.sq
edf <- mbb$edf
pval <- mbb$s.table[,"p-value"]
## list of labels
labs <- list(bquote(italic(R)^2 == .(format(r2, dig = 2))),
     bquote("edf" == .(format(edf, dig = 1))),
     bquote(italic(p) == .(format(pval, dig = 3))))
## line heights (R^2 is taller than 'edf' ...)
yoff_vec <- c(0, 1.9, 1.5)
for (i in 1:3) {
  plotrix::corner.label(labs[[i]], x=1, y=1, 
      yoff = cumsum(yoff_vec)[i]*strheight("m"))
}

I tried setting this up with atop() from ?plotmath, but it messes around with the font sizes ... I don't know if there's an easier/better way to stack lines or create a multi-line (>2) label with the plotmath machinery ...

Take 2: back to using legend() (don't need plotrix, and corner.label doesn't do "top", "middle". Use cex to set type size. Everything up to yoff_vec is the same.

lfun <- function(pos, lab, cex = 1) {
  legend(pos, legend = lab, seg.len = 0, x.intersp = 0, bty="n",
         cex = cex)
}
lfun("topleft", labs[[1]])
lfun("top",     labs[[2]])
lfun("topright", labs[[3]])

Need GAM coefficients orderly within graphs in a grid

1 Answers