Prefiro fazer análises de potência além do básico por simulação. Com pacotes pré-selecionados, nunca tenho certeza de que suposições estão sendo feitas.
A simulação da energia é bastante direta (e acessível) usando R.
decida como você deve parecer seus dados e como analisá-los
escreva uma função ou conjunto de expressões que simule os dados para um determinado relacionamento e tamanho da amostra e faça a análise (é preferível uma função, você pode transformar o tamanho e os parâmetros da amostra em argumentos para facilitar a tentativa de valores diferentes). A função ou código deve retornar o valor p ou outra estatística de teste.
use o replicate função para executar o código acima várias vezes (geralmente começo cerca de 100 vezes para ter uma idéia do tempo que leva e obter a área geral correta, depois suba para 1.000 e às vezes 10.000 ou 100.000 para o valores finais que vou usar). A proporção de vezes que você rejeitou a hipótese nula é o poder.
refaça o acima para outro conjunto de condições.
Aqui está um exemplo simples com regressão ordinal:
@gung: seu comentário faz sentido, você se importaria de adicionar seus códigos para que menos pessoas em R também pudessem se beneficiar? graças
1
Estou analisando isso novamente e tenho algumas perguntas: 1) Por que x é uniforme em 1:10? 2) Como você o generalizaria para mais de uma variável independente?
Peter Flom - Restabelece Monica
1
@ PeterFlom, x tinha que ser alguma coisa, então eu escolhi (arbitrariamente) tê-lo uniforme entre 0 e 10, também poderia ser normal, gama, etc. O melhor seria escolher algo que seja semelhante ao que esperamos do real x variáveis para se parecer. Para usar mais de uma variável preditora, gere-as independentemente (ou a partir de uma cópula normal, multivariada, etc.) e inclua-as todas na parte eta1, por exemplo eta1 <- beta0 + beta1*x1 + beta2*x2 + beta3*x3.
Greg Snow
1
replicatemeanα=0.05
3
Eu acrescentaria outra coisa à resposta de Snow (e isso se aplica a qualquer análise de potência por simulação) - preste atenção se você está procurando um teste de 1 ou 2 caudas. Programas populares como G * Power assumem o padrão de teste unicaudal e, se você estiver tentando verificar se suas simulações são compatíveis com eles (sempre é uma boa idéia quando você está aprendendo a fazer isso), verifique isso primeiro.
Para fazer do Snow um teste unilateral, eu adicionaria um parâmetro chamado "cauda" às entradas da função e colocaria algo assim na própria função:
#two-tail test
if (tail==2) fit$stats[5]
#one-tail test
if (tail==1){
if (fit$coefficients[5]>0) {
fit$stats[5]/2
} else 1
A versão unicaudal basicamente verifica se o coeficiente é positivo e, em seguida, corta o valor p pela metade.
Além do excelente exemplo de Snow, acredito que você também pode fazer uma simulação de energia reamostrando de um conjunto de dados existente que tenha seu efeito. Não é exatamente um bootstrap, já que você não está amostrando com substituição o mesmo n , mas a mesma idéia.
Então, aqui está um exemplo: fiz uma pequena auto-experiência que resultou em uma estimativa pontual positiva, mas porque era pequena, não foi quase estatisticamente significativa na regressão logística ordinal. Com esse ponto e orçamentos, como um grande n eu precisaria? Para vários n possíveis , muitas vezes gerei um conjunto de dados e executei a regressão logística ordinal e vi como o valor de p era pequeno :
library(boot)
library(rms)
npt <- read.csv("http://www.gwern.net/docs/nootropics/2013-gwern-noopept.csv")
newNoopeptPower <- function(dt, indices) {
d <- dt[sample(nrow(dt), n, replace=TRUE), ] # new dataset, possibly larger than the original
lmodel <- lrm(MP ~ Noopept + Magtein, data = d)
return(anova(lmodel)[7])
}
alpha <- 0.05
for (n in seq(from = 300, to = 600, by = 30)) {
bs <- boot(data=npt, statistic=newNoopeptPower, R=10000, parallel="multicore", ncpus=4)
print(c(n, sum(bs$t<=alpha)/length(bs$t)))
}
Nesse caso, em n = 600, a potência era de 32%. Não é muito encorajador.
(Se minha abordagem de simulação estiver errada, por favor, alguém me diga. Estou publicando alguns artigos médicos que discutem a simulação de potência para o planejamento de ensaios clínicos, mas não tenho certeza sobre minha implementação precisa.)
Ainda não tenho certeza de como a simulação deve parecer com mais (especificamente, três) variáveis independentes. Entendo que eu deveria 'incluí-los todos na parte eta1, por exemplo, eta1 <- beta0 + beta1 * x1 + beta2 * x2 + beta3 * x3' '(como mencionado acima). Mas eu não sei como ajustar o restante dos parâmetros na função. Alguém poderia me ajudar com isso?
eta1 <- beta0 + beta1*x1 + beta2*x2 + beta3*x3
.replicate
mean
Eu acrescentaria outra coisa à resposta de Snow (e isso se aplica a qualquer análise de potência por simulação) - preste atenção se você está procurando um teste de 1 ou 2 caudas. Programas populares como G * Power assumem o padrão de teste unicaudal e, se você estiver tentando verificar se suas simulações são compatíveis com eles (sempre é uma boa idéia quando você está aprendendo a fazer isso), verifique isso primeiro.
Para fazer do Snow um teste unilateral, eu adicionaria um parâmetro chamado "cauda" às entradas da função e colocaria algo assim na própria função:
A versão unicaudal basicamente verifica se o coeficiente é positivo e, em seguida, corta o valor p pela metade.
fonte
Além do excelente exemplo de Snow, acredito que você também pode fazer uma simulação de energia reamostrando de um conjunto de dados existente que tenha seu efeito. Não é exatamente um bootstrap, já que você não está amostrando com substituição o mesmo n , mas a mesma idéia.
Então, aqui está um exemplo: fiz uma pequena auto-experiência que resultou em uma estimativa pontual positiva, mas porque era pequena, não foi quase estatisticamente significativa na regressão logística ordinal. Com esse ponto e orçamentos, como um grande n eu precisaria? Para vários n possíveis , muitas vezes gerei um conjunto de dados e executei a regressão logística ordinal e vi como o valor de p era pequeno :
Com a saída (para mim):
Nesse caso, em n = 600, a potência era de 32%. Não é muito encorajador.
(Se minha abordagem de simulação estiver errada, por favor, alguém me diga. Estou publicando alguns artigos médicos que discutem a simulação de potência para o planejamento de ensaios clínicos, mas não tenho certeza sobre minha implementação precisa.)
fonte
Referindo-se à primeira simulação (sugerida por Snow; /stats//a/22410/231675 ):
Ainda não tenho certeza de como a simulação deve parecer com mais (especificamente, três) variáveis independentes. Entendo que eu deveria 'incluí-los todos na parte eta1, por exemplo, eta1 <- beta0 + beta1 * x1 + beta2 * x2 + beta3 * x3' '(como mencionado acima). Mas eu não sei como ajustar o restante dos parâmetros na função. Alguém poderia me ajudar com isso?
fonte