Sum of margins in a categorical variable with 6 levels

[Felipe Quartucci]

Hi,
I am using Stata 16.1.
I am studying the effect of different methods of land restoration on soil organic carbon sequestration. The soil organic carbon was evaluated in 6 depths (6 soil layers). I also included the soil clay content as a continuous covariate. So I have 2 independent categorical variables (method and depth - I am treating depth as categorical here), and 1 continuous independent variable (clay). Soil organic carbon (SOC) is my dependent variable.
I ran the regression and estimated the margins. The results of the margins were for each layer. But in my case, I need also the margins for the total soil organic carbon, which is the sum of all layers. So, my question is: is there a way to have a result of the margins for the sum of the layers (1+2+3+4+5+6)?

Commands used

reg SOC i.method##i.depth c.clay
contrast i.method##i.depth
margins i.method#i.depth

PS: I have already tried to include a "7th" soil depth to my data, with the sum of soil organic carbon of all other layers. But the results do not match. For example, if I sum all layers in a method, the result is not the same as the depth 7 (which was included).

[Rich Goldstein]

it is not clear to me what you want - if you want an overall test (such as one would get from anova), then you need to add "overall" to the margins contrast statement; if instead you want to save the estimates and do something else with them, Stata saves them in r(b) and the variance-covariance matrix is r(V); see

h margins

for what is saved and available to you; see

help margins_contrast
or
h contrast

[Clyde Schechter]

I'm not sure I understand what you are asking. And I don't have a clear picture of how your data is organized, but on the assumption that you have several different plots of land that have been subjected to the various treatments under study, and for each plot you have 6 observations, one for each depth layer. I'm also assuming that each plot received the same treatment at all depths. If you want to analyze total SOC for each plot, I don't think you can get that from the analyses you have shown. Rather, you would need to do something like

collapse (sum) SOC (first) method, by(plot)
regress SOC i.method

The above code ignores the clay covariate. If clay does not depend on depth and is constant within plot, then adding clay to the -collapse- command, just after method, and adding c.clay to the -regress- command will do that. But if clay varies by depth, then you need to decide what is the most sensible way to represent it in the overall analysis and modify the code accordingly.

I am not surprised that your attempted solution adding a "7th" soil depth with the sum of organic carbon to the data set failed--I would not have expected it to produce sensible results at all. Regression just doesn't work that way. I have other concerns about even the 6-levels analysis in that you are implicitly relying on the assumption that the SOC in the different depths of the same plot are independently distributed, which, even to somebody with no knowledge in this area, seems wildly implausible. I would think that a multi-level (in the statistical sense, not in the sense of depths) model would be required to capture those dependencies. In any event, to thoroughly comment on this would require more information about your data design and a view of some example data, as well as more scientific context.

Be advised that I know absolutely nothing about this subject matter and I am just trying to advise on a possible way to use Stata here. If what I am suggesting makes no sense from the perspective of soil science, then certainly disregard what I am saying.

Added: Crossed with #2.

[Felipe Quartucci]

Thank you for your attention. I tried both suggestions, but still did not get what I want. Below is how my data look like.
If I run the regression reg SOC i.method##i.depth c.clay, I will have an estimate for each layer. But I also want an estimate for the whole soil profile, which is the sum of all layers.

method	rep	depth	SOC	clay
1	1	1	26,0	407
1	1	2	31,9	435
1	1	3	73,1	461
1	1	4	40,8	484
1	1	5	28,1	484
1	1	6	32,3	506
1	2	1	36,3	204
1	2	2	23,4	359
1	2	3	42,5	410
1	2	4	22,2	360
1	2	5	11,7	153
1	2	6	9,3	153
1	3	1	15,6	151
1	3	2	15,8	151
1	3	3	40,8	175
1	3	4	21,7	201
1	3	5	9,3	176
1	3	6	10,7	150
2	1	1	49,0	382
2	1	2	44,2	458
2	1	3	63,1	482
2	1	4	39,5	488
2	1	5	34,4	513
2	1	6	52,0	487
2	2	1	40,3	202
2	2	2	24,8	434
2	2	3	49,4	483
2	2	4	31,2	460
2	2	5	17,1	462
2	2	6	23,3	538
2	3	1	27,7	177
2	3	2	19,7	228
2	3	3	40,2	226
2	3	4	26,0	252
2	3	5	12,9	254
2	3	6	12,6	100
3	1	1	50,4	487
3	1	2	41,1	561
3	1	3	82,9	585
3	1	4	58,3	585
3	1	5	40,8	562
3	1	6	46,4	597
3	2	1	92,3	307
3	2	2	40,3	567
3	2	3	62,8	568

[Clyde Schechter]

But I also want an estimate for the whole soil profile, which is the sum of all layers.

I don't understand what this means. Do you want the sum of the coefficients across all layers? You can get that with:

lincom _cons + 2.depth + 3.depth + 4.depth + 5.depth + 6.depth

If you mean that you want to combine the data from the individual layers into some overall SOC outcome, how would you do that--do you add them up? If so, after the regression you have already done for layer-by-layer results run:

collapse (sum) SOC, by(rep method)
regress SOC i.method

Or perhaps you mean something altogether different?