Custom-define contrast - mix between dummy and Helmert coding

Question

I'm trying to use custom-defined contrasts. They are sort of a combination of traditional dummy coding and the last contrast produced by reversed Helmert coding. In short, I want to compare each of the levels A, B, C with Contr (similar to dummy coding where Contr would be the reference category. Then I want a contrast comparing (A+B+C) to Contr (similar to reversed Helmert coding). I think I can define these by hand, but when trying to do it programatically in R I seem to produce something more similar to effects coding.

set.seed(123)

custom_coding <- matrix(0, ncol=4, nrow=4)
colnames(custom_coding) <- c("A-Contr","B-Contr","C-Contr",
                             "All-Contr")
rownames(custom_coding) <- c("Contr","A","B","C")

> custom_coding
      A-Contr B-Contr C-Contr All-Contr
Contr       0       0       0         0
A           0       0       0         0
B           0       0       0         0
C           0       0       0         0

custom_coding["Contr",] <- -1

# A vs Cont
custom_coding["A", grep(x=colnames(custom_coding), pattern="A-Contr")] <- 1
# B vs Contr
custom_coding["B", grep(x=colnames(custom_coding), pattern="B-Contr")] <- 1
# C vs Contr
custom_coding["C", grep(x=colnames(custom_coding), pattern="C-Contr")] <- 1
# All vs Contr
custom_coding[c("A","B","C"), grep(x=colnames(custom_coding), pattern="All-Contr")] <- 1

> custom_coding
      A-Contr B-Contr C-Contr All-Contr
Contr      -1      -1      -1        -1
A           1       0       0         1
B           0       1       0         1
C           0       0       1         1

What I would like to achieve my goal, I think is this:

contrast_matrix <- 
matrix(rbind(c(-1,1,0,0),
          c(-1,0,1,0),
          c(-1,0,0,1),
          c(-3/3,1/3,1/3,1/3)),ncol=4,
          dimnames=list(c("A-Contr","B-Contr","C-Contr","All-Contr")))

          [,1]      [,2]      [,3]      [,4]
A-Contr     -1 1.0000000 0.0000000 0.0000000
B-Contr     -1 0.0000000 1.0000000 0.0000000
C-Contr     -1 0.0000000 0.0000000 1.0000000
All-Contr   -1 0.3333333 0.3333333 0.3333333

Is there a way R can produce this automatically?

What I'm doing below, seems to result in effects coding which is not what I want.

dat <- bind_cols(data.frame(y=c(rnorm(n=10, mean=20, 3),
                                rnorm(n=10, mean=22, 3),
                                rnorm(n=10, mean=18, 3),
                                rnorm(n=10, mean=0, 5))),
      data.frame(type = c(rep(LETTERS[1:3], each=10), 
      rep("Contr", 10)))) %>% 
      mutate(type=relevel(factor(type), ref="Contr"))

contrasts(dat$type) <- solve(t(custom_coding))

> dat$type
 [1] A     A     A     A     A     A     A     A     A     A     B     B     B     B     B     B     B     B     B     B     C     C     C     C     C     C     C     C     C     C    
[31] Contr Contr Contr Contr Contr Contr Contr Contr Contr Contr
attr(,"contrasts")
      A-Contr B-Contr C-Contr
Contr    -0.5    -0.5    -0.5
A         0.5    -0.5    -0.5
B        -0.5     0.5    -0.5
C        -0.5    -0.5     0.5
Levels: Contr A B C

Even when using the previously defined contrast_matrix that I think is correct, the contrast for All-Contr drops, why?

##UPDATE I clearly missed the k-1 constraint.

> contrasts(dat$type) <- t(contrast_matrix)
> dat$type
 [1] A     A     A     A     A     A     A     A     A     A     B     B     B     B     B     B     B     B     B     B     C     C     C     C     C     C     C     C     C     C    
[31] Contr Contr Contr Contr Contr Contr Contr Contr Contr Contr
attr(,"contrasts")
      A-Contr B-Contr C-Contr
Contr      -1      -1      -1
A           1       0       0
B           0       1       0
C           0       0       1
Levels: Contr A B C

Answer 1

You are trying to define two contrasts with one matrix. To fit a model you have to choose one contrast, so that the model understands how to interpret the coefficients and estimate them.

Once you fit the model you can look at any number of comparisons between (functions of) the coefficients. For lm, glm, gls and geepack::geese linear models, you don't have to define the contrasts by hand; use the contrast package instead.

# Fit the model with the default treatment contrast. (Any contrast will do.)
fit <- lm(y ~ type, data = dat)

library("contrast")

# Example 1: Compare the control to each of the treatments.
print(
  contrast(
    fit,
    list(type = c("A", "B", "C")),
    list(type = "Contr")
  ),
  X = TRUE
)
#> lm model parameter contrast
#> 
#>  Contrast     S.E.    Lower    Upper     t df Pr(>|t|)
#>  18.61365 1.277391 16.02299 21.20432 14.57 36        0
#>  21.01564 1.277391 18.42497 23.60631 16.45 36        0
#>  15.11610 1.277391 12.52543 17.70677 11.83 36        0
#> 
#> Contrast coefficients:
#>  (Intercept) typeA typeB typeC
#>            0     1     0     0
#>            0     0     1     0
#>            0     0     0     1

# Example 2: Compare the control to the average of the treatments.
print(
  contrast(fit,
    list(type = c("A", "B", "C")),
    list(type = "Contr"),
    type = "average"
  ),
  X = TRUE
)
#> lm model parameter contrast
#> 
#>   Contrast     S.E.    Lower    Upper    t df Pr(>|t|)
#> 1 18.24847 1.042985 16.13319 20.36374 17.5 36        0
#> 
#> Contrast coefficients:
#>   (Intercept)     typeA     typeB     typeC
#> 1           0 0.3333333 0.3333333 0.3333333

# Example 3: Compare the control to a weighted average of the treatments.
print(
  contrast(fit,
    list(type = c("A", "B", "C")),
    list(type = "Contr"),
    weights = c(1 / 2, 1 / 4, 1 / 4),
    type = "average"
  ),
  X = TRUE
)
#> lm model parameter contrast
#> 
#>   Contrast     S.E.    Lower    Upper     t df Pr(>|t|)
#> 1 18.33976 1.059157 16.19169 20.48783 17.32 36        0
#> 
#> Contrast coefficients:
#>   (Intercept) typeA typeB typeC
#> 1           0   0.5  0.25  0.25

As the examples show, a contrast is a linear combination of the factor levels; it can be represented as a row vector. If you can extract the coefficients and their covariance matrix from the fitted model, you can find the contrast and its standard error by hand.

# Compute by hand the contrast that compares the control 
# to the average of the treatments.
X <- matrix(c(0, 1 / 3, 1 / 3, 1 / 3), nrow = 1)

# Find the value of the contrast
X %*% fit$coefficients
#>          [,1]
#> [1,] 18.24847
# ... and its standard error.
sqrt(X %*% vcov(fit) %*% t(X))
#>          [,1]
#> [1,] 1.042985