1
$\begingroup$

Here is an example in Zuur's book about GAM. The data file (Bailey fisheries data) and R code (Chapter2.R) can be found in his website. Below I describe what I confuse.

We want to exam whether the fish density/depth relationship is different in the two time periods. The two periods are: (1) 1979 ~ 1989; (2) 1997 ~ 2002. The GAM model is:

$Dens_{i} = \alpha + f(Depth_{i}) + \beta_{1}*Period_{i} + \epsilon_{i}$

where $Dens_{i}$ is fish density at site i, $Depth_{i}$ is the depth at site i.

The R code is:

library(mgcv)
DF <- structure(list(Dens = c(0.002070281, 0.003519799, 0.000980515, 
                              0.008039216, 0.005933375, 0.021800502, 0.009190601, 0.017616309, 
                              0.013991031, 0.013992395, 0.011129268, 0.013258145, 0.009683579, 
                              0.005310263, 0.001345291, 0.005411765, 0.001559342, 0.006637335, 
                              0.002589566, 0.008381643, 0.004223759, 0.002945137, 0.030923947, 
                              0.012547051, 0.008346952, 0.020997029, 0.014393939, 0.005137845, 
                              0.008197407, 0.004422836, 0.013695652, 0.004651163, 0.008304069, 
                              0.012608696, 0.001974744, 0.005313283, 0.006771755, 0.00468636, 
                              0.006508102, 0.020765027, 0.004004329, 0.005639414, 0.001340876, 
                              0.013311511, 0.002530903, 0.000871569, 0.023303903, 0.015236388, 
                              0.005058453, 0.012090396, 0.016701516, 0.004582651, 0.004542216, 
                              0.00056341, 0.00706278, 0.003383459, 0.008430493, 0.015, 0.004134199, 
                              0.007039963, 0.009156863, 0.000890508, 0.011035073, 0.003716609, 
                              0.00420805, 0.00325188, 0.004105058, 0.001012188, 0.009571608, 
                              9.04e-05, 0.000772639, 0.005269058, 0.008375157, 0.006525376, 
                              0.006148055, 0.002279489, 0.002211067, 0.003488906, 0.001357637, 
                              0.004353476, 0.002941587, 0.000588309, 0.003596177, 0.002487396, 
                              0.002501878, 0.00648289, 0.002549552, 0.002111706, 0.001508485, 
                              0.001508485, 0.003536285, 0.004593528, 0.003191489, 0.003441477, 
                              0.004508694, 0.004280327, 0.003322034, 0.00089252, 0.002617612, 
                              0.002376043, 0.002599356, 0.001702727, 0.000778733, 0.001800232, 
                              0.002631079, 0.000933435, 0.002495996, 0.001994553, 0.000691389, 
                              0.001020023, 0.000574192, 0.000959945, 6.15e-05, 0.001041124, 
                              0.000490196, 0.000935673, 0.000562972, 0.000505133, 0.000221903, 
                              0.001089461, 0.000332144, 0.000106167, 0.000772765, 0.000643146, 
                              0.000452546, 0.000647577, 0.00023071, 1.48e-05, 0.000575954, 
                              0.000508383, 0.00046763, 5.01e-05, 0.000288582, 0.00017508, 0.000250965, 
                              0.000140376, 0.000381388, 0.000998099, 7.78e-05, 6.79e-05, 2.65e-05, 
                              8.62e-05, 0.000155916, 8.16e-05, 0.000264026, 0.000273973, 0.00021658
), MeanDepth = c(804L, 808L, 809L, 848L, 853L, 960L, 977L, 982L, 
                 985L, 986L, 1000L, 1016L, 1017L, 1024L, 1025L, 1027L, 1036L, 
                 1053L, 1065L, 1100L, 1110L, 1144L, 1200L, 1205L, 1212L, 1217L, 
                 1257L, 1263L, 1265L, 1272L, 1284L, 1300L, 1312L, 1360L, 1368L, 
                 1379L, 1390L, 1430L, 1440L, 1448L, 1452L, 1462L, 1500L, 1506L, 
                 1519L, 1521L, 1524L, 1527L, 1533L, 1541L, 1587L, 1600L, 1642L, 
                 1648L, 1650L, 1677L, 1700L, 1720L, 1742L, 1789L, 1790L, 1845L, 
                 1867L, 1872L, 1884L, 1900L, 1909L, 1927L, 1932L, 1944L, 1946L, 
                 1975L, 1987L, 1993L, 2015L, 2016L, 2058L, 2078L, 2114L, 2122L, 
                 2172L, 2204L, 2230L, 2287L, 2292L, 2397L, 2410L, 2420L, 2434L, 
                 2450L, 2462L, 2477L, 2486L, 2487L, 2500L, 2504L, 2530L, 2572L, 
                 2645L, 2670L, 2737L, 2855L, 2875L, 2970L, 3048L, 3066L, 3089L, 
                 3099L, 3126L, 3138L, 3485L, 3639L, 3677L, 3753L, 3810L, 3900L, 
                 3985L, 3995L, 4016L, 4020L, 4056L, 4062L, 4073L, 4118L, 4127L, 
                 4222L, 4226L, 4242L, 4245L, 4292L, 4298L, 4320L, 4375L, 4510L, 
                 4565L, 4765L, 4787L, 4800L, 4810L, 4811L, 4812L, 4812L, 4840L, 
                 4842L, 4846L, 4846L, 4865L), Period = c(1L, 2L, 2L, 1L, 2L, 1L, 
                                                         1L, 1L, 1L, 1L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 1L, 2L, 2L, 
                                                         2L, 1L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 1L, 1L, 2L, 1L, 1L, 2L, 
                                                         1L, 1L, 1L, 2L, 2L, 1L, 2L, 2L, 1L, 1L, 1L, 2L, 1L, 1L, 2L, 2L, 
                                                         1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 1L, 1L, 1L, 
                                                         2L, 1L, 1L, 1L, 1L, 2L, 2L, 1L, 2L, 1L, 1L, 2L, 1L, 2L, 2L, 1L, 
                                                         2L, 1L, 2L, 2L, 2L, 1L, 1L, 1L, 2L, 2L, 1L, 2L, 1L, 1L, 1L, 1L, 
                                                         1L, 1L, 2L, 1L, 1L, 2L, 2L, 2L, 1L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 
                                                         1L, 2L, 1L, 1L, 1L, 1L, 2L, 1L, 2L, 1L, 1L, 1L, 2L, 1L, 1L, 1L, 
                                                         1L, 2L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L)), class = "data.frame", row.names = c(1L, 
                                                                                                                                                   2L, 3L, 4L, 5L, 6L, 7L, 8L, 9L, 10L, 11L, 12L, 13L, 14L, 15L, 
                                                                                                                                                   16L, 17L, 18L, 19L, 20L, 21L, 22L, 23L, 24L, 25L, 26L, 27L, 28L, 
                                                                                                                                                   29L, 30L, 31L, 32L, 33L, 34L, 35L, 36L, 37L, 38L, 39L, 40L, 41L, 
                                                                                                                                                   42L, 43L, 44L, 45L, 46L, 47L, 48L, 49L, 50L, 51L, 52L, 53L, 54L, 
                                                                                                                                                   55L, 56L, 57L, 58L, 59L, 60L, 61L, 62L, 63L, 64L, 65L, 66L, 67L, 
                                                                                                                                                   68L, 69L, 70L, 71L, 72L, 73L, 74L, 75L, 76L, 77L, 78L, 79L, 80L, 
                                                                                                                                                   81L, 82L, 83L, 84L, 85L, 86L, 87L, 88L, 89L, 90L, 91L, 92L, 93L, 
                                                                                                                                                   94L, 95L, 96L, 97L, 98L, 99L, 100L, 101L, 102L, 103L, 104L, 105L, 
                                                                                                                                                   106L, 107L, 108L, 109L, 110L, 111L, 112L, 113L, 114L, 115L, 116L, 
                                                                                                                                                   117L, 118L, 120L, 121L, 122L, 123L, 124L, 125L, 126L, 127L, 128L, 
                                                                                                                                                   129L, 130L, 131L, 132L, 133L, 134L, 135L, 136L, 137L, 138L, 139L, 
                                                                                                                                                   140L, 141L, 142L, 143L, 144L, 145L, 146L, 147L, 148L))
DF$fPeriod <- factor(DF$Period)
M1 <- gam(Dens ~ s(MeanDepth) + fPeriod, data = DF)
summary(M1)

The result shows that

The estimated density value for Period 2 is 0.0021 lower than that of Period 1.

This is also seen in "2" in the word "fPeriod2" in the summary output.

Parametric coefficients:
              Estimate Std. Error t value Pr(>|t|)    
(Intercept)  0.0054572  0.0004291  12.716  < 2e-16 ***
fPeriod2    -0.0021859  0.0007433  -2.941  0.00383 ** 

My question is: how do we determine which factor level is the "benchmark"? Because if the output is written as "fPeriod" instead of "fPeriod2", I will not know where the number "-0.0021859" should be added to (Period1 or Period2).

$\endgroup$
  • $\begingroup$ Are you sure you don’t have a typo? Because if you didn’t type fPeriod2 somewhere in your code, I don’t see why R would give you a coefficient for it. $\endgroup$ – Joe Aug 11 '19 at 3:16
  • 1
    $\begingroup$ @Joe It seems not a typo. Maybe the "mgcv" package (Version 1.8.24) is really humanizing. $\endgroup$ – T X Aug 12 '19 at 8:30
  • $\begingroup$ @Joe It's because "The name of this dummy variable is obtained by pasting together the word Period and the non-reference factor level 2: Period2". See Isabella Ghement's answer (should be "fPeriod" by the way). My last comment hasn't seen the answer $\endgroup$ – T X Aug 12 '19 at 8:38
5
$\begingroup$

When R includes a factor in a regression-type model such as yours, it automatically converts that factor to a set of dummy variables.

R accomplishes this by setting aside the first level of the factor - which will be treated as the reference level - and introducing dummy variables for comparing the mean value of the response variable for all other levels against this reference level (with the comparison being adjusted for the effects of the other predictors in the model).

When creating these dummy variables, R assigns names to them so that it can keep track of their effects when reporting the model summary. The name of each dummy variable is simply obtained by pasting together the name of the factor and the name of the non-reference factor level corresponding to that dummy variable.

In your example, Period has two levels: 1 and 2. R sets asside the first level (1) as a reference and compares the second level (2) against this reference level via a dummy variable. The name of this dummy variable is obtained by pasting together the word Period and the non-reference factor level 2: Period2.

What does R consider to be the first factor level? You can find this out using the levels() command:

levels(DF$Period)

The first level listed by R for the Period factor in the R console is what R will by default consider as the reference level of Period.

You can change the reference level of a factor using the relevel() command:

DF$Period <- relevel(DF$Period, ref = "2")

So you don't have to do any typing yourself for R to create names for dummy variables it creates behind the scenes to encode the effect of a factor - it just automatically does the naming for you, as explained above.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.