How to deal with bias caused by data collected with different methods on a GLM? I have two datasets from different years (2020 and 2021) that I needed to merge because my sample was too small. I'm investigating which morphological features (sex, age, weight, fat reserves and muscle mass) from a migratory bird species affect the number of days they're staying in a determined area before leaving for migration (they have GPS tags on them).
The problem is that the 2020 data was collected by someone else and 2021 was collected by me. While all the variables are equal and measured right, the person before me captured and measured birds in different months while I opted to capture and measure all birds on the same month (as soon as birds start to aggregate in the study area). Because it wasn't my initial plan to merge the datasets, this caused problems.
I'm using negative binomial distribution because my response variable is count data (number of days) and on R the model is glm.nb = days ~ age + sex + fat + muscle + weight. The summary of this model presents me significant values for muscle, fat and age, which match my initial hypothesis but as soon as a put year on the mix (a variable I feel I can't ignore), every other variable loses it's significance and year is the only one significant. Nagelkerke's R² is bigger on the model with the year variable and the AIC is lower, but I'm sure the way data was collected is causing this.
Is there something I can do to eliminate this effect of "year" or is it just a lost cause? I have the exact dates of each capture so this may be better than year? I'm new here so sorry if more information is needed or I'm not being clear.
 A: 
I have the exact dates of each capture so this may be better than year?

This might provide a solution if you break it down by month as a predictor, with the month in which you collected as the reference level. You presumably have overlap with the other investigator for the calendar month during which you collected data. You would still include year as a predictor, with the year of the other investigator as a reference level. Then the year coefficient would be for the (covariate-adjusted) difference between the two estimates for that month, while coefficients for the other month categories would be the differences from the reference month during the year that the other investigator did the collections.
Second, consider whether including weight as a predictor is really adding anything. With my limited knowledge of ornithology, it seems that would mostly be the sum of fat plus muscle and perhaps correlated with age, so it would be highly correlated to other predictors. Even for predictors that are actually important, individual coefficient estimates for correlated predictors can seem "insignificant" with high variances but counterbalancing coefficient-estimate covariances.
Third, think about whether you need to model days as count data. It seems that it might better be modeled as a continuous outcome.
Finally, if the number of days is only the time they stayed in the region after they were captured and tagged, then you have censored data. In that case you only know the minimum number of days that they stayed; the actual time could have been longer. That would require a regression analysis that takes censoring into account, like survival analysis. If you or the other investigator did any tagging, that might help explain the discrepancy associated with year.
