2012-02-03 1 views
0

직교 대비가 c1과 c2 인 3 그룹 범주 형 예측 변수 x를 정규 분포 변수 y라고 가정 해 봅시다. 나는 R에 프로그램을 만들려고합니다. x, c1, c2가 주어지면, c1과 c2는 사용자가 지정한 효과 크기 r1과 r2를 갖도록 y를 만듭니다. 예를 들어주어진 모집단 효과 크기로 ANOVA 형 디자인에서 무작위 추출

,의 X가, C1은 C2는, R1, R2는 같이 생성 된 것을 가정 해 봅시다 다음

x <- factor(rep(c(1, 2, 3), 100)) 
contrasts(x) <- matrix(c(0, -.5, .5, -2/3, 1/3, 1/3), 
    nrow = 3, ncol = 2, dimnames = list(c("1", "2", "3"), c("c1", "c2"))) 

contrasts(x) 
    c1   c2 
1 0.0 -0.6666667 
2 -0.5 0.3333333 
3 0.5 0.3333333 

r1 <- .09 
r2 <- 0 

나는 Y의 분산이 차지하도록 Y를 만들 수있는 프로그램을하고 싶습니다 c1은 r1 (.09)이고 y의 분산은 c2가 차지하는 r2 (0)과 같습니다.

아무도 내가이 문제를 어떻게 알 수 있습니까? 나는 내가 rnorm 함수를 사용해야한다는 것을 알고 있지만, 어떤 인구 집단이/sds rnorm이 그것의 샘플링을 할 때 사용해야 하는지를 고집합니다. 당신의 도움에 미리

감사합니다,

패트릭 내 동료의 일부 관대 한 조언

답변

0

씨는, 지금, 그룹의 지정된 번호, 대조의 집합 주어진 시뮬레이션 데이터를 생성 한 기능을 가지고 회귀 계수의 세트 셀당 지정된 N 및 소정 그룹 내 분산

sim.factor <- function(levels, contr, beta, perCell, errorVar){ 
    # Build design matrix X 
    X <- cbind(rep(1,levels*perCell), kronecker(contr, rep(1,perCell))) 
    # Generate y 
    y <- X %*% beta + rnorm(levels*perCell, sd=sqrt(errorVar)) 
    # Build and return data frame 
    dat <- cbind.data.frame(y, X[,-1]) 
    names(dat)[-1] <- colnames(contr) 
    return(dat) 
} 

는 I는 회귀 계수, 셀당 N, 기수 세트 소정의 집합 함수를 작성 직교 반대로 STS 원하는 델타 R 관심 콘트라스트^2, 필요한 집단 내 분산 리턴로는 다음 몇 가지 테스트를 수행 한 후

ws.var <- function(levels, contr, beta, perCell, dc){ 
    # Build design matrix X 
    X <- cbind(rep(1,levels), contr) 
    # Generate the expected means 
    means <- X %*% beta 
    # Find the sum of squares due to each contrast 
    var <- (t(means) %*% contr)^2/apply(contr^2/perCell, 2, sum) 
    # Calculate the within-conditions sum of squares 
    wvar <- var[1]/dc - sum(var) 
    # Convert the sum of squares to variance 
    errorVar <- wvar/(3 * (perCell - 1)) 
    return(errorVar) 
} 

를 함수는 콘트라스트 C1의 원하는 델타 R^2를 생성하는 것 .

contr <- contr.helmert(3) 
colnames(contr) <- c("c1","c2") 
beta <- c(0, 1, 0) 
perCell <- 50 
levels = 3 
dc <- .08 
N <- 1000 

# Calculate the error variance 
errorVar <- ws.var(levels, contr, beta, perCell, dc) 

# To store delta R^2 values 
d1 <- vector("numeric", length = N) 

# Use the functions 
for(i in 1:N) 
{ 
    d <- sim.factor(levels=3, 
        contr=contr, 
        beta=beta, 
        perCell=perCell, 
        errorVar=errorVar) 
    d1[i] <- lm.sumSquares(lm(y ~ c1 + c2, data = d))[1, 2] # From the lmSupport package 
} 

m <- round(mean(d1), digits = 3) 

bmp("Testing simulation functions.bmp") 
hist(d1, xlab = "Percentage of variance due to c1", main = "") 
text(.18, 180, labels = paste("Mean =", m)) 
dev.off() 

패트릭

관련 문제