Your code works fine (apart from the typo) if you (naïvely) assume a normal distribution (i.e. a t-distribution with infinite degrees of freedom).
all.equal(qnorm(1 - .05/2), qt(1 - .05/2, df=Inf))
# [1] TRUE
`colnames<-`(t(apply(d, 1, function(x)
x[1] + x[2]*(qt(1 - .05/2, df=Inf)*c(-1, 1)))), paste0(c(2.5, 97.5), "%"))
# 2.5% 97.5%
# Gender -1.08798019 1.06798019
# Age -0.06879892 0.04879892
# Job satisfaction -0.27679712 0.03679712
However, you have n=112 observations, m=3 coefficients and k=1 constant, thus n - m - k degrees of freedom. Hence using the t-distribution with 108 degrees of freedom might be the better choice.
(DOF <- 112 - 3 - 1)
# [1] 108
qt(1 - .05/2, df=DOF)
# [1] 1.982173
`colnames<-`(t(apply(d, 1, function(x)
x[1] + x[2]*(qt(1 - .05/2, df=DOF)*c(-1, 1)))), paste0(c(2.5, 97.5), "%"))
# 2.5% 97.5%
# Gender -1.1001954 1.08019542
# Age -0.0694652 0.04946520
# Job satisfaction -0.2785739 0.03857388
For a complete summary you could add t-statistics and p-values
signif(cbind(d, t=d[,1]/d[,2], p=2*pt(-abs(d[,1]/d[,2]), df=DOF),
`colnames<-`(
t(apply(d, 1, function(x)
x[1] + x[2]*(qt(1 - .05/2, df=DOF)*c(-1, 1)))),
paste0(c(2.5, 97.5), "%"))),
2)
# B.hat. se t p 2.5% 97.5%
# Gender -0.01 0.55 -0.018 0.99 -1.100 1.100
# Age -0.01 0.03 -0.330 0.74 -0.069 0.049
# Job satisfaction -0.12 0.08 -1.500 0.14 -0.280 0.039
Data:
d <- structure(list(B.hat. = c(-0.01, -0.01, -0.12), se = c(0.55,
0.03, 0.08)), class = "data.frame", row.names = c("Gender", "Age",
"Job satisfaction"))