Estimando porcentagens como variável dependente na regressão

Eu tenho as porcentagens de classificação dos alunos em 38 exames como a variável dependente no meu estudo. Uma porcentagem de classificação é calculada por (classificação / número de alunos em um exame). Essa variável dependente tem distribuição quase uniforme e eu quero estimar os efeitos de algumas variáveis ​​na variável dependente.

Qual abordagem de regressão eu uso?

Se você estiver trabalhando com a Stata, consulte o seguinte exemplo: http://www.ats.ucla.edu/stat/stata/faq/proportion.htm

Aqui está uma citação desta página da web:

"Como se faz regressão quando a variável dependente é uma proporção?

Os dados da proporção têm valores que variam entre zero e um. Naturalmente, seria bom que os valores previstos também caíssem entre zero e um. Uma maneira de conseguir isso é usar um modelo linear generalizado (glm) com um link de logit e a família binomial. Incluiremos a opção robusta no modelo glm para obter erros padrão robustos que serão particularmente úteis se tivermos especificado incorretamente a família de distribuição ".

Sinopse

Os resultados da regressão podem ter algum valor limitado quando cuidadosamente interpretados. Formas inevitáveis ​​de variação farão com que as estimativas do coeficiente diminuam substancialmente para zero. É necessário um modelo melhor que lide com a variação de maneira mais apropriada.

(Um modelo de máxima verossimilhança pode ser construído, mas pode ser impraticável devido ao cálculo necessário, o que envolve a avaliação numérica de integrais multidimensionais. O número de dimensões é igual ao número de alunos matriculados nas aulas.)

Introdução

Como narrativa para fundamentar nossa intuição, imagine que esses 38 exames foram realizados em 38 cursos separados durante um semestre em uma pequena escola com matrícula de 200 estudantes universitários. Em uma situação realista, esses alunos terão diferentes habilidades e experiências. Como medidas substitutas dessas habilidades e experiências, podemos obter, digamos, notas nos testes de matemática e verbais do SAT e ano na faculdade (1 a 4).

Normalmente, os alunos se matriculam em cursos de acordo com suas habilidades e interesses. Os calouros fazem cursos introdutórios e os cursos introdutórios são preenchidos principalmente por calouros. Os alunos do ensino médio e os calouros e alunos do segundo ano do ensino médio fazem os cursos avançados e de pós-graduação. Essa seleção estratifica parcialmente os alunos, de modo que as habilidades inatas dos alunos de qualquer classe são tipicamente mais homogêneas do que a disseminação de habilidades por toda a escola.

Assim, os alunos mais capazes podem encontrar-se pontuando perto da parte inferior das aulas difíceis e avançadas nas quais se matriculam, enquanto os alunos menos capazes podem pontuar perto do topo das aulas introdutórias fáceis que recebem. Isso pode confundir uma tentativa direta de relacionar as classificações dos exames diretamente aos atributos dos alunos e das turmas.

Análise

Indexe os alunos com e deixe que os atributos do aluno i sejam dados pelo vetor x i . Indexe as classes com j e deixe que os atributos da classe j sejam dados pelo vetor z j . O conjunto de alunos matriculados na classe j é A $i$$i$$\mathbf{x}_i$$j$$j$$\mathbf{z}_j$$j$$A_j$ .

Assumir a "força" de cada aluno é uma função de seus atributos mais algum valor aleatório, que pode muito bem ter média zero:$s_i$

Modelamos o exame na classe adicionando valores aleatórios independentes à força de cada aluno matriculado na classe e convertendo-os em classificações. Daí, se estudante i é inscrito na classe j , a sua posição relativa r i , j é determinada pela sua posição na matriz ordenada de valores$j$$i$$j$$r_{i,j}$

Esta posição é dividido por um mais do que o total de matrículas classe para dar a variável dependente, a ordem percentual:$r_{i,j}$

Afirmo que os resultados da regressão dependem (um pouco) dos tamanhos e da estrutura dos valores aleatórios (não observados) e δ i , j . $\varepsilon_i$$\delta_{i,j}$ Os resultados também dependem precisamente de como os alunos estão matriculados nas aulas. Isso deve ser intuitivamente óbvio, mas o que não é tão óbvio - e parece difícil de analisar teoricamente - é como e quanto os valores não observados e as estruturas de classe afetam a regressão.

Simulação

Sem muito esforço , podemos simular essa situação para criar e analisar alguns dados de amostra. Uma vantagem da simulação é que ela pode incorporar o verdadeiro pontos fortes dos alunos, que na realidade não são observáveis. Outra é que podemos variar os tamanhos típicos dos valores não observados, bem como as atribuições da classe. Isso fornece uma "caixa de areia" para avaliar os métodos analíticos propostos, como a regressão.

Para começar, vamos definir o gerador de números aleatórios para resultados reproduzíveis e especificar o tamanho do problema. Eu uso Rporque está disponível para qualquer pessoa.

set.seed(17)
n.pop <- 200      # Number of students
n.classes <- 38   # Number of classes
courseload <- 4.5 # Expected number of classes per student

Para fornecer realismo, crie n.classesclasses de dificuldades variadas em duas escalas (matemática e verbal, com correlação negativa), conduzidas em diferentes níveis acadêmicos (variando de 1 = introdutório a 7 = pesquisa) e com facilidade variável. (Em uma aula "fácil", as diferenças entre as quantidades de aprendizado dos alunos podem ser grandes e / ou o exame pode oferecer pouca discriminação entre os alunos. Isso é modelado por termos aleatórios que, para a classe j, tendem a ser grandes Os resultados do exame serão quase imprevisíveis a partir dos dados de força do aluno. Quando a turma não é "fácil", esses termos aleatórios são insignificantes pequenos e a força do aluno pode determinar perfeitamente a classificação do exame.)$\delta_{i,j}$$j$

classes <- data.frame(cbind(
  math <- runif(n.classes), 
  rbeta(n.classes, shape1=(verbal <- (1-math)*5), shape2=5-verbal),
  runif(n.classes, min=0, max=7),
  rgamma(n.classes, 10, 10)))
rm(math, verbal)
colnames(classes) <- c("math.dif", "verbal.dif", "level", "ease")
classes <- classes[order(classes$math.dif + classes$verbal.dif + classes$level), ]
row.names(classes) <- 1:n.classes
plot(classes, main="Classes")

Os alunos são distribuídos entre os quatro anos e recebem valores aleatórios de seus atributos. Não há correlações entre nenhum destes atributos:

students <- data.frame(cbind(
  as.factor(ceiling(runif(n.pop, max=4))),
  sapply(rnorm(n.pop, mean=60, sd=10), function(x) 10*median(c(20, 80, floor(x)))),
  sapply(rnorm(n.pop, mean=55, sd=10), function(x) 10*median(c(00, 80, floor(x)))),
  rnorm(n.pop)
  ))
colnames(students) <- c("year", "math", "verbal", "ability")
plot(students, main="Students")

O modelo é que cada aluno tem uma "força" inerente determinada em parte por seus atributos e em parte por sua "habilidade", que é o valor . Os coeficientes de força , que determinam a força em termos de outros atributos, são o que a análise de dados subsequente procurará estimar. Se você quiser jogar com esta simulação, faça-o alterando . A seguir, é apresentado um conjunto interessante e realista de coeficientes, refletindo o aprendizado contínuo dos alunos durante a faculdade (com uma grande quantidade entre os anos 2 e 3); onde 100 pontos em cada parte do SAT valem cerca de um ano de escola; e onde cerca de metade da variação é devida aos valores de "habilidade" não capturados pelas notas do SAT ou pelo ano na escola.$\varepsilon_i$betabeta

beta <- list(year.1=0, year.2=1, year.3=3, year.4=4, math=1/100, verbal=1/100, ability=2, sigma=0.01)
students$strength <- (students$year==1)*beta$year.1 + 
  (students$year==2)*beta$year.2 +
  (students$year==3)*beta$year.3 +
  (students$year==4)*beta$year.4 +
  students$math*beta$math + 
  students$verbal*beta$verbal + 
  students$ability*beta$ability
students <- students[order(students$strength), ]
row.names(students) <- 1:n.pop

(Ter em mente que students$ability é inobservável: é um desvio aparentemente aleatório entre a força prevista dos outros atributos observáveis ​​e a força real nos exames. Para remover esse efeito aleatório, defina beta$abilitycomo zero. beta$sigmaMultiplicará os easevalores: é basicamente o desvio padrão do relação à gama de forças dos alunos em um determinado curso Valores em torno de 0,01 a $\delta_{i,j}$$.01$$.2$ ou mais me parecem razoáveis.)

Deixe os alunos escolherem cursos que correspondam às suas habilidades. Uma vez que eles façam isso, podemos calcular os tamanhos das classes e esconder esses com o classesdataframe para uso posterior. O valor de spreadna assignments <-...linha determina o quanto os alunos são divididos em classes por habilidade. Um valor próximo de combina essencialmente os alunos mais fracos com os cursos mais fáceis. Um valor próximo ao número de aulas espalha os alunos um pouco mais. Valores muito maiores do que esses começam a ficar irreais, porque tendem a colocar os alunos mais fracos nos cursos mais difíceis.$0$

pick.classes <- function(i, k, spread) {
  # i is student strength rank
  # k is number to pick
  p <- pmin(0.05, diff(pbeta(0:n.classes/n.classes, i/spread, (1+n.pop-i)/spread)))
  sample(1:n.classes, k, prob=p)
}
students$n.classes <- floor(1/2 + 2 * rbeta(n.pop,10,10) * courseload)
assignments <- lapply(1:n.pop, function(i) pick.classes(i, students$n.classes[i], spread=1))
enrolment <- function(k) length(seq(1, n.pop)[sapply(assignments, function(x) !is.na(match(k, x)))])
classes$size <- sapply(1:n.classes, enrolment)
classes$variation <- by(data, data$Class, function(x) diff(range(x$strength)))

(Como um exemplo do que esta etapa foi realizada, veja a figura mais abaixo.)

Agora aplique o modelo: as habilidades dos alunos em cada turma são variadas de forma independente - mais para exames fáceis, menos para exames difíceis (discriminatórios) - para determinar suas notas. Eles são resumidos como classificações e "brincadeiras", que são porcentagens de classificação. As brincadeiras para uma classe de alunos variam de 1 / ( n + 1 ) a n / ( n $n$$1/(n+1)$$n/(n+1)$$1/(n+1)$$0$$1$

exam.do <- function(k) {
  s <- seq(1, n.pop)[sapply(assignments, function(x) !is.na(match(k, x)))]
  e <- classes$ease[k]
  rv <- cbind(rep(k, length(s)), s, order(rnorm(length(s), students$strength[s], sd=e*beta$sigma*classes$variation[k])))
  rv <- cbind(rv, rv[,3] / (length(s)+1))
  dimnames(rv) <- list(NULL, c("Class", "Student", "Rank", "Prank"))
  rv
}
data.raw <- do.call(rbind, sapply(1:n.classes, exam.do))

A esses dados brutos, anexamos os atributos de aluno e classe para criar um conjunto de dados adequado para análise:

data <- merge(data.raw, classes, by.x="Class", by.y="row.names")
data <- merge(data, students, by.x="Student", by.y="row.names")

Vamos nos orientar inspecionando uma amostra aleatória dos dados:

> data[sort(sample(1:dim(data)[1], 5)),]

Row Student Class Rank Prank math.dif verbal.dif  level  ease Size year math verbal ability strength n.classes
118      28     1   22 0.957  0.77997   6.95e-02 0.0523 1.032   22    2  590    380   0.576     16.9         4
248      55     5   24 0.889  0.96838   1.32e-07 0.5217 0.956   26    3  460    520  -2.163     19.0         5
278      62     6   22 0.917  0.15505   9.54e-01 0.4112 0.497   23    2  640    510  -0.673     19.7         4
400      89    10   16 0.800  0.00227   1.00e+00 1.3880 0.579   19    1  800    350   0.598     21.6         5
806     182    35   18 0.692  0.88116   5.44e-02 6.1747 0.800   25    4  610    580   0.776     30.7         4

O registro 118, por exemplo, diz que o aluno nº 28 se matriculou na classe nº 1 e ficou em 22º (de baixo) no exame para uma classificação percentual de 0,957. O nível geral de dificuldade desta classe foi de 0,0523 (muito fácil). Um total de 22 alunos foram matriculados. Este estudante está no segundo ano (segundo ano) com 590 de matemática, 380 notas verbais no SAT. Sua força acadêmica inerente é de 16,9. Eles estavam matriculados em quatro classes na época.

Este conjunto de dados comporta com a descrição na pergunta. Por exemplo, as classificações percentuais são quase uniformes (como devem ser para qualquer conjunto de dados completo, porque as classificações percentuais para uma única classe têm uma distribuição uniforme e discreta).

Lembre-se, em virtude dos coeficientes em beta, este modelo assumiu uma forte conexão entre as pontuações dos exames e as variáveis ​​mostradas neste conjunto de dados. Mas o que mostra a regressão? Vamos regredir a logística da classificação percentual em relação a todas as características observáveis ​​dos alunos que possam estar relacionadas às suas habilidades, bem como aos indicadores de dificuldade da classe:

logistic <- function(p) log(p / (1-p))
fit <- lm(logistic(Prank) ~ as.factor(year) + math + verbal + level, data=data)
summary(fit)

Coefficients:
                  Estimate Std. Error t value Pr(>|t|)    
(Intercept)      -2.577788   0.421579   -6.11  1.5e-09 ***
as.factor(year)2  0.467846   0.150670    3.11   0.0020 ** 
as.factor(year)3  0.984671   0.164614    5.98  3.2e-09 ***
as.factor(year)4  1.109897   0.171704    6.46  1.7e-10 ***
math              0.002599   0.000538    4.83  1.6e-06 ***
verbal            0.002130   0.000514    4.14  3.8e-05 ***
level            -0.208495   0.036365   -5.73  1.4e-08 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

Residual standard error: 1.48 on 883 degrees of freedom
Multiple R-squared: 0.0661, Adjusted R-squared: 0.0598 
F-statistic: 10.4 on 6 and 883 DF,  p-value: 3.51e-11 

As plotagens de diagnóstico ( plot(fit)) parecem fastásticas: os resíduos são homocedásticos e lindamente normais (embora um pouco curtos, o que não é problema); sem discrepâncias; e nenhuma influência desagradável em qualquer observação.

$3.5$$(-9, 1.6, 3.4, 3.9, 0.009, 0.007, -0.7)$$(*, 1, 3, 4, 0.010, 0.010, *)$$*$ significa um coeficiente que não foi especificado explicitamente).

level$0$$7$level

(A propósito, o uso das classificações percentuais não transformadas na regressão não altera qualitativamente os resultados relatados abaixo.)

spread$1$$38$

Coefficients:
                  Estimate Std. Error t value Pr(>|t|)    
(Intercept)      -4.902006   0.349924  -14.01  < 2e-16 ***
as.factor(year)2  0.605444   0.130355    4.64  3.9e-06 ***
as.factor(year)3  1.707590   0.134649   12.68  < 2e-16 ***
as.factor(year)4  1.926272   0.136595   14.10  < 2e-16 ***
math              0.004667   0.000448   10.41  < 2e-16 ***
verbal            0.004019   0.000434    9.25  < 2e-16 ***
level            -0.299475   0.026415  -11.34  < 2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

Residual standard error: 1.3 on 883 degrees of freedom
Multiple R-squared: 0.282,  Adjusted R-squared: 0.277 
F-statistic: 57.9 on 6 and 883 DF,  p-value: <2e-16

( Neste gráfico de dispersão de tarefas da turma, spreaddefinido como$38$spread1

Desta vez, o R-quadrado é muito melhorado (embora ainda não seja ótimo). No entanto, todos os coeficientes aumentaram de 20 a 100%. Esta tabela os compara com algumas simulações adicionais:

Simulation Intercept Year.2 Year.3 Year.4 Math Verbal Level R^2
Beta               *    1.0    3.0    4.0 .010   .010     *   *
Spread=1        -2.6    0.5    1.0    1.1 .003   .002 -0.21  7%
Spread=38       -4.9    0.6    1.7    1.9 .005   .004 -0.30 25%
Ability=1       -8.3    0.9    2.6    3.3 .008   .008 -0.63 58%
No error       -11.2    1.1    3.3    4.4 .011   .011 -0.09 88%

Guardando spread$38$ability$2$$1$abilitysigma$0$$\varepsilon_i$$\delta_{i,j}$ completamente, dá um elevado ao quadrado R e produz estimativas perto com os valores correctos. (Vale ressaltar que o coeficiente levelentão diminui em uma ordem de magnitude.)

Essa análise rápida mostra que a regressão, pelo menos como realizada aqui, confunde formas inevitáveis ​​de variação com os coeficientes. Além disso, os coeficientes também dependem (até certo ponto) de como os alunos são distribuídos entre as aulas. Isso pode ser parcialmente acomodado pela inclusão de atributos de classe entre as variáveis ​​independentes na regressão, como feito aqui, mas mesmo assim o efeito da distribuição dos alunos não desaparece.

Qualquer falta de previsibilidade do verdadeiro desempenho do aluno, e qualquer variação no aprendizado do aluno e no desempenho real nos exames, aparentemente fazem com que as estimativas do coeficiente diminuam para zero. Eles parecem fazê-lo uniformemente, sugerindo que os coeficientes relativos ainda possam ser significativos.

A medida que o usuário13203 propõe pode ser considerada como uma pontuação de desempenho inferior limitado, quanto menor, melhor o desempenho: $y_{ij}$ i-ésimo desempenho do aluno no j-ésimo exame.

Usando uma transformação linearizante de logit em que $\mu_{ij}$ pode depender das características observáveis ​​dos alunos ou dos exames:

$\ln(y_{ij}/(1-y_{ij})) = \mu_{ij} + e_{ij} + v_i$

as habilidades não observadas do aluno são modeladas através do componente aleatório $v_i$ enquanto $e_{ij}$modela outros não observáveis ​​não sistemáticos. A correlação entre respostas (exames) pode ser corrigida afirmando uma estrutura geral de covariância para$e_{ij}$. Por que não uma estrutura de variação de branco (ou sanduíche / robusta)? Além disso, algumas das correlações de respostas podem ser explicadas dentro do$\mu_{ij}$ (dependência condicional).

(Esta é apenas uma ideia da minha experiência tendenciosa, comentários e críticas são mais que bem-vindos.)

É provável que habilidades não observáveis ​​sejam correlacionadas com atributos observáveis ​​de alunos ou exames dentro $\mu_{ij}$. Essa premissa faz deste modelo um ER com componentes de erro correlacionados, que podem ser estimados por ML ou por um estimador de dois estágios: primeiro estágio: uma transformação interna (ou analógica) que elimina$v_i$. Segunda etapa: OLS no modelo transformado.

Você pode tentar a regressão logística. A transformação de logit$\ln(\frac{p}{1-p})$ espalhará sua variável de resposta pela linha real para que você não obtenha porcentagens de classificação previstas absurdas como -3% ou + 110%.

Um modelo perfeito nesse caso mapeará as entradas (quaisquer que sejam suas covariáveis) para as saídas (a classificação do aluno na classe). Outra maneira de pensar nisso é mapeando primeiro as pontuações e, em seguida, mapeando essas pontuações para a classificação. Vou ignorar o erro por enquanto.

pontuação do teste: $y = \sum \beta x$

classificação: $r = R(y)$

No qual $R$é a função de classificação. O problema é que$R$é uma função não linear que depende inteiramente dos próprios dados. Se assumirmos que temos uma quantidade infinita de dados, saberemos a distribuição completa de$y$e $R(y)$é essencialmente a função de densidade cumulativa. Ele informa qual porcentagem de pessoas obteve uma pontuação pior que você no teste, a área à esquerda da sua pontuação.

Isso parece ser bastante semelhante à forma funcional do modelo linear generalizado. Penso que é por isso que a abordagem de regressão logística foi proposta por Mike Anderson. Se as pontuações dos seus exames forem distribuídas logisticamente, a função de link a ser usada será o logit (sua inversa é a função de densidade cumulativa com a qual nos preocupamos). Da mesma forma, se as pontuações fossem normalmente distribuídas, a função probit seria a função de link.

Para sua regressão, a única maneira de estimar as classificações é dizer "dado que meus dados são distribuídos como X, esse ponto está no 34º percentil". Caso contrário, como você sabe o que significa um aumento de dois pontos na pontuação do teste em termos de classificação? A ressalva é que você precisa estimar essa distribuição para escolher sua função de link (certas formas funcionais tornarão sua vida muito mais fácil). Além disso, esse modelo não diz "você foi o sexto melhor de uma classe de 38", e sim "se as pontuações dos testes fossem distribuídas como pensamos que são, sua pontuação o colocaria no percentil 15".