[R-meta] Question about meta-regression with multiple treatment studies

[Mensching, Andre]

Dear all experts in meta-analaysis,

I'm currently working on a meta-analysis in agricultural science. Unfortunately, I am very uncertain about the method to be used. After recently submitting a paper, I received major revisions and the advice to use the metafor package and in particular the rma.mv function. So far I only used linear mixed models and the lmer function from the lme4 package in R and the “inverse-variance” method.

I think it would be best to briefly outline the objective and describe the data structure:
I want to do a meta regression on the subject of metabolic diseases in dairy cows. The target variable is the daily mean pH value (= continues outcome) over a period of e.g. 48 h in the forehead system. The published studies are mostly based on replicated Latin Square experimental designs. For example 8 cows are split over 4 different treatments (= 2 x 2 factorial design with different diets, a "control group" is not necessarily present) within 1 of 4 periods. In the different periods the cows receive different treatments, thus each of the 8 cows is exposed to all treatments.  
The studies provide LSmean estimates per treatment for the mean pH as well as the standard error of these estimates (not the standard deviation!). In addition, information on treatments in form of the feed composition and feed analyses are provided. This information as well as milk constituents, for example, are to be examined as moderator variables.

Here is an example for the data structure in the long-format (simulated random numbers, DMI (=Dry Matter Intake) and NDF (Neutral Detergent Fiber) are moderators for the diet, additionally the milk fat content):

study pHmean pHmean_SE n_cows       DMI      NDF      fat
        1           6.09           0.05                     8    23.75    18.88     4.01
        1           6.08           0.05                     8    23.60    17.71     3.93
        1           5.84           0.05                     8    23.16    18.51     3.94
        1           5.96           0.05                     8    25.18    17.93     3.91
        2           5.76           0.09                   16    25.98    19.53     3.91
        2           6.00           0.09                   16    23.29    19.30     4.01
        2           6.13           0.09                   16    27.72    20.57     3.97
        2           5.92           0.09                   16    22.17    18.61     3.99
        3           5.98           0.09                     8    26.50    19.54     3.91
        3           6.00           0.09                     8    24.76    20.74     3.95
        3           6.12           0.09                     8    23.35    20.58     4.02
        3           5.86           0.09                     8    27.76    19.73     4.08 
       …               ...                …                    ...           …           ...        …  

I think that it is necessary to work with multiple mixed-effects models and to consider the raw means as effect size of each treatment. 

As already mentioned, so far I used lmer() [study coded as factor!) :

LMM <- lmer(pHmean ~ DMI + NDF + fat  + (1|study), 
	weights =  1/pHmean_SE, 
      	REML = TRUE, 
      	data = Data)

I also tried the lme function from the package nlme:

LMM2 <- lme(pHmean ~ DMI + NDF + fat ,                               
  	random = ~ 1 | study,
 	weights = varFixed(~pHmean_SE),
  	control=lmeControl(sigma = 1),
  	method="REML",
  	data=Data)

My first approach with the rma.mv function is:

res <- rma.mv(yi = pHmean, 
	V = pHmean_SE,  
     	mods = ~ DMI + NDF + fat  ,
                random = ~ 1 | study,
  	method="REML",
                struct="UN",  
	data=Data)

To the theme is also already a statement:
http://www.metafor-project.org/doku.php/tips:rma_vs_lm_lme_lmer 

I can confirm that if nlme and lmer return the same result, if the control argument is not set to lmeControl(sigma = 1). If I set the argument control=lmeControl(sigma = 1), I get similar results to those of the rma.mv function (but not the same!).

During my literature research on studies with meta-analyses using a comparable data structure, I only found evaluations using either lmer (R) or proc mixed (SAS) with weights = 1/SEM  or  1/SEM^2. 

For example:
Roman-Garcia, Y., R.R. White, and J.L. Firkins. 2016. Meta-analysis of postruminal microbial nitrogen flows in dairy cattle. I. Derivation of equations. J. Dairy Sci. 99:7918–7931. http://dx.doi.org/10.3168/jds.2015-10661.
Santos, J.E.P., I.J. Lean, H. Golder, and E. Block. 2019. Meta-analysis of the effects of prepartum dietary cation-anion difference on performance and health of dairy cows. J. Dairy Sci. 102:2134–2154. http://dx.doi.org/10.3168/jds.2018-14628.
White, R.R., M.B. Hall, J.L. Firkins, and P.J. Kononoff. 2017. Physically adjusted neutral detergent fiber system for lactating dairy cow rations. I: Deriving equations that identify factors that influence effectiveness of fiber. J. Dairy Sci. 100:9551–9568. http://dx.doi.org/10.3168/jds.2017-12765.

To summarize the problems again: 
I actually used a method that others have used and even recently was published.
Nevertheless I was advised to use the metafor package with the rma.mv function, because it should be more appropriate. 
And yes, the results differ of course. So I do not know what is right and what is not. 

What do you think about that?

I would be very happy to hear from someone. Thank you very much in advance!

Best regards,
André