Aplicação do Trabalho de Conclusão de Curso intitulado como “Evasão e conclusão no curso de bacharelado em estatística: uma aplicação da análise de sobrevivência na presença de riscos competitivos”

O presente arquivo trata-se de um recorte das análises do TCC , cujo objetivo geral é analisar a evasão e a conclusão no curso de Bacharelado em Estatística no Brasil, identificando variáveis internas e externas que influenciem tais fenômenos. Além disso, são adotados alguns objetivos específicos tais quais: traçar o perfil dos alunos que estão mais propensos a experimentar os eventos de interesse; identificar e quantificar as variáveis que aceleram ou retardam o tempo até a evasão ou a conclusão e predizer as funções de incidência acumulada para perfis selecionados de indivíduos.

• Dados Utilizados

Foram utilizados os dados do Censo da Educação Superior (Censup) que, de acordo com o Instituto Nacional de Estudos e Pesquisas Educacionais Anísio Teixeira (Inep), é o instrumento de pesquisa mais completo do Brasil acerca da Educação Superior desde 1995.

É possível acessar os microdados do Censup no sítio do Inep e as informações podem ser obtidas via download, em formato ASCII.

Utilizou-se o Censup do ano de 2010 referente aos alunos que ingressaram no curso de Bacharelado em Estatística no mesmo ano. Para acompanhar a situação desses alunos na IES, foram considerados os anos de 2010 a 2017.

As variáveis código do aluno e código da situação do aluno foram utilizadas para conectar as bases de dados do Censup do ano de 2010 com os anos de 2011 até 2017, com o intuito de traçar a trajetória destes alunos no curso ao longo dos anos. Isto foi feito utilizando a variável que indica o código da situação do aluno, que neste estudo foi recategorizada.

Considerou-se o aluno ativo, aquele que estava cursando ou com a matrícula trancada; evadido, aquele que havia se desvinculado do curso ou transferido para outro curso na mesma IES e, por fim, concluso, aquele aluno que havia se formado.

• Software R

Para o cálculo da Função de incidência acumulada, pode-se utilizar a função survifit da biblioteca survival, que estima a incidência acumulada para os eventos de interesse, utilizando o argumento etype e, sendo de interesse, pode-se incorporar covariáveis. Ou alternativamente, pode-se utilizar a biblioteca cmprsk e a função cuminc que também é possivel calcular as FIA’s.

Foi utilizada a função crprep da biblioteca mstate para calcular os pesos que cada indivíduo assume à medida que ocorrem eventos competitivos, é repetido o mesmo procedimento para cada um dos eventos de interesse, assumindo os demais como competitivos. Com isso, é criado um novo banco de dados e esse que é utilizado na modelagem.

Após o procedimento realizado para o cálculo dos pesos, com a função coxph estima-se os efeitos das covariáveis, indicando o evento de interesse no argumento failcode.

#bibliotecas

#carregando o banco de dados já com os devidos filtros de interesse

Banco de dados considerando o recorte e filtros mencionados acima.

setwd("C:\\Users\\Adriane\\Documents\\7º período\\TCC")
banco<-read.table("banco_final.txt", sep="", header = TRUE)
attach(banco)
names(banco)
##  [1] "nu"             "SITUACAO2010"   "SITUACAO2011"   "SITUACAO2012"  
##  [5] "SITUACAO2013"   "SITUACAO2014"   "SITUACAO2015"   "SITUACAO2016"  
##  [9] "SITUACAO2017"   "SITUACAO2010_1" "SITUACAO2011_1" "SITUACAO2012_1"
## [13] "SITUACAO2013_1" "SITUACAO2014_1" "SITUACAO2015_1" "SITUACAO2016_1"
## [17] "SITUACAO2017_1" "Tempo"          "Censura"        "Evento"        
## [21] "Adm"            "acad"           "raca"           "sexo"          
## [25] "reserva"        "apoio"          "turno"          "ativ"          
## [29] "vest"           "idade"

#Função de incidência acumulada (FIA) - geral

Cálculo da FIA geral, isto é, sem estratificação por alguma covariável de interesse.

#geral
layout(1)
ci1<-cuminc(Tempo, Evento, rho=0, cencode=0)
plot(ci1, lty=1,lwd=2, col=1:2, xlab="tempo(anos)", ylab="probabilidade ",ylim=c(0,0.9),
     main='', cex.main=1,
     curvlab = c("conclusão", "evasão"))

As curvas da FIA referentes aos dois eventos de interesse indicam que um aluno está mais propenso a vir a experimentar primeiro a evasão do que a conclusão. Nos primeiros anos é esperado que a FIA para a conclusão seja próxima de zero, dado que o curso de bacharelado em estatística tem duração típica de 4 anos.

Ao final de 8 anos, a probabilidade de os alunos terem concluído foi de aproximadamente 0,27, enquanto para a evasão, esta probabilidade já era superior nos primeiros anos

#FIA’s considerando covariaveis

Agora analisaremos os comportamentos das FIA’s de acordo com 8 covariáveis de interesse: sexo, idade, turno, participação em atividades complementares, apoio social, forma de ingresso, uso de reserva de vaga e tipo de organização acadêmica, a fim de traçar o perfil dos alunos mais propensos a evadir e concluir.

#por covariaveis

#sexo
sexoisi<-rep(NA,length(banco$sexo))
sexoisi[sexo=='Masculino']<-'masculino'
sexoisi[sexo=='Feminino']<-'feminino'
ci_sexisi<-cuminc(Tempo, Evento, sexoisi, rho=0, cencode=0)
plot( ci_sexisi,col=c(1,1,2,2), xlab="tempo(anos)",
      main='', cex.lab=.9,
      cex=.6, lty=c(1,2,1,2),lwd=c(2,2,2,2), ylab="probabilidade" )
title("Sexo")

#idade
idaisi<-rep(NA,length(banco$idade))
idaisi[min(idade)<=idade&idade<quantile(idade,0.5)]<-'idade<mediana'
idaisi[idade>=quantile(idade,0.5)]<-'idade>mediana'
ci_idaisi<-cuminc(Tempo, Evento, idaisi, rho=0, cencode=0)
plot(ci_idaisi,col=c(1,1,2,2), xlab="tempo(anos)",
     main='', cex.lab=.9,
     cex=.6, lty=c(1,2,1,2),lwd=c(2,2,2,2), ylab="probabilidade" )
title("Idade")

#turno
turnoisi<-rep(NA,length(banco$turno))
turnoisi[banco$turno==3]<-'noturno'
turnoisi[banco$turno<=2]<-'não noturno'
turnoisi[banco$turno>3]<-'não noturno'
ci_turnoisi<-cuminc(Tempo, Evento, turnoisi, rho=0, cencode=0)
plot(ci_turnoisi,col=c(1,1,2,2), xlab="tempo(anos)",
     main='', cex.lab=.9,
     cex=.6, lty=c(1,2,1,2),lwd=c(2,2,2,2), ylab="probabilidade" )
title("Turno")

#atividade complementar
atvisi<-rep(NA,length(banco$ativ))
atvisi[banco$ativ==1]<-'engajado em atividades'
atvisi[banco$ativ==0]<-'não engajado em atividades'
ci_atisi<-cuminc(Tempo, Evento, atvisi, rho=0, cencode=0)
plot( ci_atisi,col=c(1,1,2,2), xlab="tempo(anos)",
      main='', cex.lab=.9,
      cex=.6, lty=c(1,2,1,2),lwd=c(2,2,2,2), ylab="probabilidade" )
title("Atividade complementar")

#reserva de vaga
reser1<-rep(NA,length(banco$reserva))
reser1[banco$reserva==1]<-'reserva'
reser1[banco$reserva==0]<-'não reserva'
ci_res<-cuminc(Tempo, Evento, reser1, rho=0, cencode=0)
plot( ci_res,col=c(1,1,2,2), xlab="Tempo (anos)",
      main='', cex.lab=.9,
      cex=.6, lty=c(1,2,1,2),lwd=c(2,2,2,2), ylab="Probabilidade" )
title("Reserva de vaga")

#apoio social
apoio1<-rep(NA,length(banco$apoio))
apoio1[banco$apoio==1]<-'apoio social'
apoio1[banco$apoio==0]<-'não apoio social'
ci_ap<-cuminc(Tempo, Evento, apoio1, rho=0, cencode=0)
plot( ci_ap,col=c(1,1,2,2), xlab="Tempo (anos)",
      main='', cex.lab=.9,
      cex=.6, lty=c(1,2,1,2),lwd=c(2,2,2,2), ylab="Probabilidade" )
title("Apoio social")

#vestibular
vest1<-rep(NA,length(banco$vest))
vest1[banco$vest==1]<-'vestibular'
vest1[banco$vest==0]<-'não vestibular'
ci_vest<-cuminc(Tempo, Evento, vest1, rho=0, cencode=0)
plot( ci_vest,col=c(1,1,2,2), xlab="Tempo (anos)",
      main='', cex.lab=.9,
      cex=.6, lty=c(1,2,1,2),lwd=c(2,2,2,2), ylab="Probabilidade" )
title("Ingresso")

#organização acadêmica
acad1<-rep(NA,length(banco$acad))
acad1[banco$acad<2]<-'Universidade'
acad1[banco$acad>=2]<-'Outro tipo'
ci_acad<-cuminc(Tempo, Evento, acad1, rho=0, cencode=0)
plot( ci_acad,col=c(1,1,2,2), xlab="Tempo (anos)",
      main='', cex.lab=.9,
      cex=.6, lty=c(1,2,1,2),lwd=c(2,2,2,2), ylab="Probabilidade" )
title("Organização acadêmica")

#Modelagem

Na etapa da modelagem, é necessário primeiro criar um novo banco de dados em que incorporamos a informação dos pesos dos índividuos a medida que ocorre um dos eventos de interesse.

#criando o banco com os pesos

#categorizando#
#dim(banco)

sexoo<-rep(NA,length(banco$sexo))
sexoo[sexo=='Masculino']<-0
sexoo[sexo=='Feminino']<-1

acad1<-rep(NA,length(banco$acad))
acad1[acad<2]<-1
acad1[acad>=2]<-0

turno1<-rep(NA,length(banco$turno))
turno1[banco$turno<3]<-0
turno1[banco$turno==3]<-c(1)
turno1[banco$turno>3]<-0

ida<-rep(NA,length(banco$idade))
ida[min(idade)<=idade&idade<quantile(idade,0.5)]<-'0'
ida[idade>=quantile(idade,0.5)]<-'1'

banco$acad1<-acad1
banco$turno1<-turno1
banco$sexo1<-sexoo
banco$ida<-ida

#criando novo banco
banco.pesos<-crprep(Tstop = "Tempo", status = "Evento", 
                    trans=c(1,2), id= "nu",cens = 0, 
                    keep=c("SITUACAO2017_1","sexo1","acad1","sexo","reserva","apoio","turno1","ativ","vest","idade", "ida"),
                    data = banco)
banco.pesos[1:5,]
##   nu Tstart Tstop status weight.cens SITUACAO2017_1 sexo1 acad1      sexo
## 1  1      0     5      1           1      ConclusÒo     0     1 Masculino
## 2  2      0     2      2           1    DesistÛncia     0     1 Masculino
## 3  2      2     8      2           1    DesistÛncia     0     1 Masculino
## 4  3      0     1      2           1    DesistÛncia     0     1 Masculino
## 5  3      1     8      2           1    DesistÛncia     0     1 Masculino
##   reserva apoio turno1 ativ vest idade ida count failcode
## 1       0     0      0    0    1    19   0     1        1
## 2       0     0      0    0    1    20   1     1        1
## 3       0     0      0    0    1    20   1     2        1
## 4       0     0      0    0    1    19   0     1        1
## 5       0     0      0    0    1    19   0     2        1

Após a seleção de variáveis, obtivemos as seguintes variaveis significativas para os evetos de interesse:

Portanto, a modelagem e a estimação dos coeficientes serão para estes grupos de covariáveis.

#estimação dos coeficientes após seleção de covariaveis

#modelo final selecionado - conclusão
con_fim<-coxph(Surv(Tstart, Tstop, SITUACAO2017_1=="ConclusÒo")~sexo1+idade+acad1+turno1+ativ,
             data=banco.pesos, weights = weight.cens, subset = failcode ==1)
summary(con_fim)
## Call:
## coxph(formula = Surv(Tstart, Tstop, SITUACAO2017_1 == "ConclusÒo") ~ 
##     sexo1 + idade + acad1 + turno1 + ativ, data = banco.pesos, 
##     weights = weight.cens, subset = failcode == 1)
## 
##   n= 2655, number of events= 445 
## 
##            coef exp(coef) se(coef)      z Pr(>|z|)    
## sexo1   0.36855   1.44564  0.09585  3.845 0.000121 ***
## idade  -0.05887   0.94283  0.01072 -5.493 3.95e-08 ***
## acad1  -0.59368   0.55229  0.15640 -3.796 0.000147 ***
## turno1 -0.47000   0.62500  0.10391 -4.523 6.09e-06 ***
## ativ    1.02334   2.78246  0.18074  5.662 1.50e-08 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
##        exp(coef) exp(-coef) lower .95 upper .95
## sexo1     1.4456     0.6917    1.1980    1.7444
## idade     0.9428     1.0606    0.9232    0.9628
## acad1     0.5523     1.8106    0.4065    0.7504
## turno1    0.6250     1.6000    0.5098    0.7662
## ativ      2.7825     0.3594    1.9525    3.9653
## 
## Concordance= 0.648  (se = 0.013 )
## Likelihood ratio test= 120.7  on 5 df,   p=<2e-16
## Wald test            = 114.4  on 5 df,   p=<2e-16
## Score (logrank) test = 120.1  on 5 df,   p=<2e-16
#modelo final selecionado - desistência
eva_fim<-coxph(Surv(Tstart, Tstop, SITUACAO2017_1=="DesistÛncia")~sexo1+idade+turno1+ativ,
             data=banco.pesos, weights = weight.cens, subset = failcode ==2)
summary(eva_fim)
## Call:
## coxph(formula = Surv(Tstart, Tstop, SITUACAO2017_1 == "DesistÛncia") ~ 
##     sexo1 + idade + turno1 + ativ, data = banco.pesos, weights = weight.cens, 
##     subset = failcode == 2)
## 
##   n= 2021, number of events= 1073 
## 
##             coef exp(coef)  se(coef)      z Pr(>|z|)    
## sexo1  -0.219303  0.803079  0.064073 -3.423 0.000620 ***
## idade   0.018789  1.018966  0.004095  4.588 4.48e-06 ***
## turno1  0.232081  1.261222  0.061945  3.747 0.000179 ***
## ativ   -0.798733  0.449899  0.206640 -3.865 0.000111 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
##        exp(coef) exp(-coef) lower .95 upper .95
## sexo1     0.8031     1.2452    0.7083    0.9105
## idade     1.0190     0.9814    1.0108    1.0272
## turno1    1.2612     0.7929    1.1170    1.4240
## ativ      0.4499     2.2227    0.3001    0.6745
## 
## Concordance= 0.575  (se = 0.009 )
## Likelihood ratio test= 73.32  on 4 df,   p=5e-15
## Wald test            = 71.09  on 4 df,   p=1e-14
## Score (logrank) test = 72.5  on 4 df,   p=7e-15

O Modelo Conclusão indica que a taxa de falha da conclusão do curso de estatística para as mulheres é de 1,45 vezes a taxa de conclusão dos homens. Além disso, a cada aumento de um ano na idade, a taxa de conclusão diminui em 6%. Os alunos matriculados em uma IES do tipo universidade têm a taxa de conclusão reduzida em 45% quando comparado com o grupo dos alunos matriculados em IES de outro tipo. Para os alunos que frequentam o curso noturno a taxa de conclusão é 37% menor em relação àqueles que frequentam os demais turnos. Por fim, o aluno estar envolvido em atividades complementares no primeiro ano da graduação aumenta a taxa de conclusão em 2,78 vezes quando comparada aos alunos que não participam de tais atividades.

O Modelo Evasão mostra que a taxa de falha da evasão do curso de Estatística para as mulheres é 20% menor que a taxa de evasão para os homens. A cada aumento de um ano na idade, a taxa de evasão aumenta em 1,8%. Para os alunos que frequentam o curso noturno a taxa de evasão é 1,26 vezes a taxa de evasão para os alunos que não frequentam o curso noturno e, por fim, os alunos que fazem atividade completar no primeiro ano da graduação possui a taxa de evasão 55% menor que a taxa dos alunos que não fazem atividade complementar.

#Predição das FIA’s para perfis selecionados

Para calcular as taxas de falha de conclusão e de evasão de alunas que estudam em diferentes turnos e com diferentes engajamentos em atividades complementares, foram consideradas duas mulheres de 25 anos de idade que estudam em uma organização acadêmica de tipo não universidade. A primeira aluna está matriculada em curso não noturno e realiza atividade complementar; a segunda aluna está matriculada em curso noturno e não realiza atividade complementar.

tt<-seq(1:8)

#conclusao

ss<- survfit(con_fim)
s0<-round(ss$surv,digits=5)
H0<- -log(s0)
b1<-con_fim$coefficients[1]
b2<-con_fim$coefficients[2]
b3<-con_fim$coefficients[3]
b4<-con_fim$coefficients[4]
b5<-con_fim$coefficients[5]
#individuo com : sexo=fem ; idade=25 ; acad=n?o universidade ; turno: n?o noturno e ativ=sim
pred3=b1*1+b2*25+b3*0+b4*0+b5*1
f3<-1-exp(-H0*exp(pred3))

#individuo com : sexo=fem ; idade=25 ; acad=n?o universidade ; turno: noturno e ativ=n?o
pred4=b1*1+b2*25+b3*0+b4*1+b5*0
f4<-1-exp(-H0*exp(pred4))

# Os dois tipos de layout para rodar

layout(1)
par(mfrow=c(1,2))

plot(tt,f3,type="s",ylim=range(c(0,1)),xlab="Tempos",ylab="Probabilidade",lty=1)
lines(tt,f4,type="s",lty=2, col=1)
legend(1,1,lty=c(1,2),c("não noturno e realiza atividade","noturno e não realiza atividade"),
       lwd=c(2,2),col=c(1,1),bty="n",cex=0.7)
title("Conclusão")

#desistencia
ss<- survfit(eva_fim)
s0<-round(ss$surv,digits=5)
H0<- -log(s0)
b1<-eva_fim$coefficients[1]
b2<-eva_fim$coefficients[2]
b3<-eva_fim$coefficients[3]
b4<-eva_fim$coefficients[4]

#individuo com : sexo=fem ; idade=25 ; turno: n?o noturno e ativ=sim
pred3=b1*1+b2*25+b3*0+b4*1
f3<-1-exp(-H0*exp(pred3)) 

#individuo com : sexo=fem ; idade=25 ; turno: noturno e ativ=n?o
pred4=b1*1+b2*25+b3*1+b4*0
f4<-1-exp(-H0*exp(pred4))

#par(mfrow=c(1,1))
plot(tt,f3,type="s",ylim=range(c(0,1)),xlab="Tempos",ylab="Probabilidade",lty=1, col=2)
lines(tt,f4,type="s",lty=2, col=2)
legend(1,1,lty=c(1,2),c("não noturno e realiza atividade","noturno e não realiza atividade"),
       lwd=c(2,2),col=c(2,2),bty="n",cex=0.7)
title("Evasão")

Vemos que o fato de frequentar curso não noturno e fazer atividade complementar acelerou a taxa de conclusão. O contrário foi verificado em relação à evasão, ou seja, o fato de estar matriculado no curso noturno e não fazer atividade complementar aumentou significativamente a taxa de evasão.

#observações finais

Vale ressaltar que a aplicação do TCC aqui exposta é apenas um recorte de uma análise maior, em que nosso intuito é ilustrar o uso das principais bibliotecas e funções utilizadas no TCC na abordagem da Análise de Sobrevivência na presença de riscos competitivos, sendo que não devem se limitar apenas a essas.

Para melhor compreensão e aplicação das funções mencionadas para análise de riscos competitivos no software R, ver Carvalho et al. (2011).