I believe in that instance they are accounting for additional covariates by first regressing them out using a mixed-effects model. Then HAllA was run on the matrices of residuals instead of the matrices of raw data.
So in your example, you’d run a mixed effects model to regress out the effect of the environment variable that you think might cause confounding.
Given how dataset-specific this modelling would need to be, HAllA doesn’t include the functionality to run the model and extract the residuals. You can see more discussion about this type of process in this thread: