m2regress {biology} | R Documentation |
~~ A concise (1-5 lines) description of what the function does. ~~
m2regress(X, Y, type = "RMA", zero = F)
X |
~~Describe X here~~ |
Y |
~~Describe Y here~~ |
type |
~~Describe type here~~ |
zero |
~~Describe zero here~~ |
~~ If necessary, more details than the description above ~~
~Describe the value returned If it is a LIST, use
comp1 |
Description of 'comp1' |
comp2 |
Description of 'comp2' |
...
....
~~further notes~~
~Make other sections like Warning with section{Warning }{....} ~
~~who you are~~
~put references to the literature/web site here ~
~~objects to See Also as help
, ~~~
##---- Should be DIRECTLY executable !! ---- ##-- ==> Define data, use random, ##-- or do help(data=index) for the standard data sets. ## The function is currently defined as function (X, Y, type = "RMA", zero = F) { IV <- X DV <- Y ssx <- sum((IV - mean(IV))^2) ssy <- sum((DV - mean(DV))^2) sxy <- sum(((IV - mean(IV))) * ((DV - mean(DV)))) r <- sxy/sqrt(ssx * ssy) r2 <- sxy/(ssx * ssy) df.DV <- length(DV) - 2 n <- length(IV) s2xy <- abs((ssy - (sxy^2)/ssx)/(length(IV) - 2)) mod <- lm(DV ~ IV) if (type == "OLS") { b <- sxy/ssx a <- mean(DV) - b * mean(IV) pp <- length(IV) sb <- sqrt((s2xy/ssx)) sa <- sqrt(s2xy * (1/(length(IV) - 0) + (mean(IV)^2)/ssx)) ci <- qt(0.975, length(IV) - 2) * sb bLower <- b - ci bUpper <- b + ci aCI <- qt(0.975, length(IV) - 2) * sa aLower <- a - aCI aUpper <- a + aCI } if (type == "MA") { b_ols <- sxy/ssx r2 <- (sxy/sqrt(ssx * ssy))^2 d <- (b_ols^2 - r2)/(r2 * b_ols) b <- (d + sqrt(d^2 + 4))/2 a <- mean(DV) - b * mean(IV) sb <- sqrt((s2xy/ssx)) sa <- sqrt(s2xy * (1/(length(IV) - 0) + (mean(IV)^2)/ssx)) B <- ((qt(0.975, n - 2)^2)) * ((1 - r2))/(n - 2) bUpper <- b * (sqrt(B + 1) + sqrt(B)) bLower <- b * (sqrt(B + 1) - sqrt(B)) aLower <- mean(DV) - (bUpper * mean(IV)) aUpper <- mean(DV) - (bLower * mean(IV)) } if (type == "rMA") { if (zero == F) { rIV <- (X - min(X))/(max(X) - min(X)) rDV <- (Y - min(Y))/(max(Y) - min(Y)) } if (zero == T) { rIV <- (X)/(max(X)) rDV <- (Y)/(max(Y)) } rssx <- sum((rIV - mean(rIV))^2) rssy <- sum((rDV - mean(rDV))^2) rsxy <- sum(((rIV - mean(rIV))) * ((rDV - mean(rDV)))) rr <- rsxy/sqrt(rssx * rssy) rr2 <- rr^2 rb_ols <- rsxy/rssx d <- (rb_ols^2 - rr2)/(rr2 * rb_ols) b <- (d + sqrt(d^2 + 4))/2 rn <- length(rIV) rs2xy <- abs((rssy - (rsxy^2)/rssx)/(length(rIV) - 2)) if (zero == F) b <- b * ((max(Y) - min(Y))/(max(X) - min(X))) if (zero == T) b <- b * ((max(Y))/(max(X))) a <- mean(DV) - b * mean(IV) sb <- sqrt((rs2xy/rssx)) sa <- sqrt(rs2xy * (1/(length(rIV) - 0) + (mean(rIV)^2)/rssx)) B <- ((qt(0.975, rn - 2)^2)) * ((1 - rr2))/(n - 2) bUpper <- b * (sqrt(B + 1) + sqrt(B)) bLower <- b * (sqrt(B + 1) - sqrt(B)) aLower <- mean(DV) - (bUpper * mean(IV)) aUpper <- mean(DV) - (bLower * mean(IV)) } if (type == "RMA") { r2 <- (sxy/sqrt(ssx * ssy))^2 b <- sqrt(ssy/ssx) if (sxy < 0) b <- b * -1 a <- mean(DV) - (b * mean(IV)) sb <- sqrt((s2xy/ssx)) sa <- sqrt(s2xy * (1/(length(IV) - 0) + (mean(IV)^2)/ssx)) n <- length(IV) ci <- qt(0.975, length(IV) - 2) * sb bLower <- b - ci bUpper <- b + ci B <- ((qt(0.975, n - 2)^2)) * ((1 - r2))/(n - 2) bUpper <- b * (sqrt(B + 1) + sqrt(B)) bLower <- b * (sqrt(B + 1) - sqrt(B)) aCI <- qt(0.975, length(IV) - 2) * sa aLower <- a - aCI aUpper <- a + aCI aLower <- mean(DV) - (bUpper * mean(IV)) aUpper <- mean(DV) - (bLower * mean(IV)) } coef <- matrix(NA, 2, 1, dimnames = list(c("(Intercept)", "X"), c("Coef"))) stats <- matrix(NA, 2, 5, dimnames = list(c("(Intercept)", "X"), c("SE", "t", "P|>t|", "2.5", "97.5"))) coef[1, 1] <- a coef[2, 1] <- b stats[1, 1] <- sa stats[2, 1] <- sb stats[1, 2] <- a/sa stats[2, 2] <- b/sb stats[1, 3] <- 2 * pt(abs(abs(a/sa)), length(DV) - 2, lower.tail = F) stats[2, 3] <- 2 * pt(abs(b/sb), length(DV) - 2, lower.tail = F) stats[1, 4] <- aLower stats[1, 5] <- aUpper stats[2, 4] <- bLower stats[2, 5] <- bUpper z <- list(coefficients = coef, stats = stats) z }