Classificando amostras

FAZ UMA CLASSIFICACAO DAS AMOSTRAS EM GRUPOS A PRIORI A PARTIR DOIS PASSOS:

1. Executa um teste de classificação para as todas as amostras que tem um categoria de grupo atribuida (leave-one-out do SCRIPT-02)
2. Duas alternativas:
   1. escolher o modelo com maior porcentagem de acerto global (svm, bayes ou lda), definir um modelo discriminante com amostras preditas corretamente e reclassificar as demais com base nesse modelo, refinando a classificação, obtendo a probabilidade de pertenecimento à classe predita. Meio circular, talvez. Primeiro script abaixo.
   2. para todas as amostras classificadas corretamente pelas n/3 análises efetuadas (svm, bayes, lda) atribui a categoria GRUPO CORE, i.e. amostras que define os grupos. Com dados dessas amostras faz um modelo (svm, bayes ou lda) e prediz a categoria para as demais amostras, classificando elas e obtendo a probabilidade de pertenecimento à classe predita.
   3. As probabilidade podem ser usadas para expressar incerteza na classificação. As amostras do grupo core tem 100% de probabilidade, pois definem o modelo em primeiro lugar.

Testa a classificação inicial (a priori) das amostras para identificar amostras para um modelo mais refinado. Requer função testa.classificacao(), disponível em SCRIPT-02 Testando grupos.

library("e1071")
dad = dadosMORFOnmds #contém os dados morfológicos e uma coluna com as classes iniciais
colgrp = "GROUPING" #nome da coluna que tem a classificacao inicial
dad[,colgrp] = as.factor(dad[,colgrp])
#testa a classificação
oteste = testa.classificacao(dad,colgrp=colgrp)  
#pega os resutados brutos com o resultados combinados de svm, bayes e lda 
brutos = oteste$BRUTOS

#ALTERNATIVA 01 - UTILIZA APENAS O TESTE COM MELHOR RESULTADO
#filtrar por apenas 1 modelo
vl = brutos$ANALISE=='SVM'
brutos = brutos[vl,]
 
#quais amostras devem fazer parte do modelo preditivo (GRUPO CORE)?
 
#definA uma probabilidade mínima de pertenecimento ao GRUPO CORE, que irá gerar o modelo preditivo com o qual a classificação de todas as amostras será feita
prob.min = 0.95 
#se for muito strito talvez não tenha n suficiente para todos os grupos
#então vamos ver primeiro como foram todos os acertos
vl = brutos$ACERTO=='acerto'
#variacao de probabilidade de pertenecimento
summary(brutos$MELHOR.PREDICAO.PROB[vl])
#numero de amostras que foram preditas corretamente (acertos) por categoria predita
table(as.vector(brutos$MELHOR.PREDICAO)[vl])
#agora com limite da probabilidade minima
vl2 = brutos$MELHOR.PREDICAO.PROB>prob.min & vl
table(as.vector(brutos$MELHOR.PREDICAO)[vl2])
 
#qual a solucao final para o modelo?
#escolhendo todas as amostras de acerto
vl = as.vector(brutos$ACERTO)=='acerto'
train.ids = as.vector(brutos$ID[vl])
test.ids = as.vector(brutos$ID[!vl]) #erros para serem reclassificados

#pega os resutados brutos 
brutos = oteste$BRUTOS
 
#pega identificadores únicos para as amostras
amostras = as.vector(unique(brutos$ID))
 
#define o criterio para selecao de amostra para GRUPO CORE (modelo) = numero de modelos, e prob.min
n.modelos =2 #precisa ter sido acerto em pelo menos 2 análises
prob.min = 0 #tanto faz a probabilidade de pertenecimento ao grupo se for acerto
 
#funcao que pegaids segundo os criterios acima
entra.no.modelo <- function(x,dd=brutos,n.modelos=1,prob.min=0) {
  #tem mais acertos entre todas as análises svm, bayes, lda?
  vl = dd$ID==x & dd$ACERTO=='acerto' & dd$MELHOR.PREDICAO.PROB>prob.min
  #print(sum(vl))
  if (sum(vl)>=n.modelos) {
    return(TRUE)
  } else { return(FALSE)}
}  
 
#que amostras entrariam no modelo caso apenas 1 modelo tenha acertado
#mas tenha acertado com prob.min>0.95
quais = sapply(amostras,entra.no.modelo,dd=brutos,n.modelos=1,prob.min=0.95)
#quantas tem por grupo
quais = amostras[quais]
table(dad[quais,colgrp])
 
#muito restrito? diminua prob.min e/ou n.modelos
 
#proporcao de amostras que irao no modelo?
length(quais)/nrow(dad)
train.ids = quais
test.ids = rownames(dad)[!rownames(dad)%in%quais]

#define um modelo com as amostras core (train.ids definidos por uma das alternativas acima)
svmm <- svm(dad[train.ids,colgrp] ~ ., dad[train.ids,], probability=TRUE)
 
#pega o dado test (pelos test.ids definidos por uma das alternativas acima
test = dad[test.ids,]
 
#classifica essas amostras segundo o modelo core
predsvm = predict(svmm,test,probability=T)
 
#qual foi a melhor predicao
preditocomo = as.vector(predsvm)
#quais sao as probabilidades de predicao de todas as categorias
predicoes = as.data.frame(attr(predsvm,"probabilities"))
 
#adiciona colunas
predicoes$ORIGINAL = original
predicoes$MELHOR.PREDICAO = preditocomo
predicoes$MELHOR.PREDICAO.PROB = NA
mytrim <- function (x) {gsub("^\\s+|\\s+$", "", x)}
for(p in 1:nrow(predicoes)) {
  pt = mytrim(preditocomo[p])
  val = predicoes[p,pt]
  predicoes$MELHOR.PREDICAO.PROB[p] = val
}
 
#anota as incertezas e junta com o modelo
 
#classe das amostras no modelo (por definicao a probabilidade é maxima)
md = data.frame(ID=train.ids,CLASSE=dad[train.ids,colgrp],GRUPO.CORE=1,PROB.CLASSE=1,stringsAsFactors = F)
 
#junta com as classes preditas
td = data.frame(ID=test.ids,CLASSE=predicoes[test.ids,"MELHOR.PREDICAO"],GRUPO.CORE=0,PROB.CLASSE=round(predicoes[test.ids,"MELHOR.PREDICAO.PROB"],2),stringsAsFactors = F)
 
#junta tudo
resultado = rbind(md,td)
 
#salva a classificacao
write.table(resultado,file="reclassificacaoSVM.csv",sep="\t",na="",row.names = F)