I am working with a dataset that consists of abundance counts of 6 microbial taxa in a lake measured weekly for 20 weeks. I also have environmental data (temperature, nutrient concentrations, conductivity etc) associated with each bacterial census. In total, i have about 25 environmental variables (IVs), though i expect only 4 or 5 of them to be strongly correlated with my DVs (species abundance).
I would like to construct a model that would relate the abundance of my six species (DVs) to environmental conditions in the lake (IVs). I expect my response variables to have unimodal responses to my dependent variables - as is most commonly the case in ecology. Typically a negative-binomial distribution is recommended for this type of species count data. I also expect there to be interactions between my IVs because the species are competing for resources. There is a great R package called mvabund which will fit GLMs to each of the species in your dataset to model their abundance. However, a major assumption is independence of the response variable. In my case, the abundance of a species should be temporally autocorrelated.
What is an appropriate model to use to relate species abundance to environmental variables in timeseries data? I am not interested in predictive power, but rather in understanding the relative importance of different environmental variables to each species' abundance.