I am trying to model the proportions of a variable with four categories through time and at four different sites. Most sites have all four categories of the variable, but some sites only have three categories. I want to model the trend in the categories across years using a hierarchical GAM that includes a shared smooth, as well as a site-specific smooth.

Consider the multinomial GAM with mgcv in R