====== Calculando as transições ======
Este script e ainda experimental para o calculo das transicoes.

<code>
### Preparando os dados o cálculo das transições ###

#Estes procedimentos foram preparados para uma planilha típica de dados de demografia de árvores, onde as linhas são os indivíduos e as colunas as medidas tomadas, sendo que repetições temporais da mesma medida representam colunas diferentes.
#O indivíduo pode ter uma das seguintes medidas: DAP(diâmetro à altura do peito), para plantas com DAP até 70mm (com paquímetro), PAP (perímetro à altura do peito - fita métrica) ou DAS (diâmetro à altura do solo), quando à altura do peito ele apresenta DAP<10mm ou tem altura inferior a 1,3m.
#As medidas em cada tempo estão identificadas por um número: 0, 1 e 2

#Preparando a planilha a ser salva como txt:
#1. Espaços em branco (principalmente nos outros tipos de medida): o R coloca NA automaticamente, mas lembrar de trocar o código "NA" (de não se aplica) no status, pois pode haver confusão. Trocar para "NE" (de não existe).
#2. Excluir observações
#3. Simplificar títulos das colunas!
#4. Não esquecer de arrumar os códigos das colunas de status

teste <- read.table("demo.gua.txt", header=T, sep="\t")

#Inserção das colunas com as áreas basais (para a somatória dos que tem caules múltiplos e dos que tem mais de um tipo de medida)

teste$ab.dap0 <- round((pi*(teste$dap0^2))/4,2)		## valores arredondados para dois dígitos
teste$ab.dap1 <- round((pi*(teste$dap1^2))/4,2)
teste$ab.das0 <- round((pi*(teste$das0^2))/4,2)		
teste$ab.das1 <- round((pi*(teste$das1^2))/4,2)
teste$ab.pap0 <- round((teste$pap0^2)/(4*pi),2)
teste$ab.pap1 <- round((teste$pap1^2)/(4*pi),2)

#Somatória das áreas basais por placa
ab.final <- aggregate(teste[,c("ab.dap0", "ab.dap1", "ab.das0", "ab.das1", "ab.pap0", "ab.pap1")], by=list(teste$placa), FUN=sum, na.rm=T)	##indexação de todas as linhas do objeto teste e as colunas especificadas
colnames(ab.final)=c("placa", "dap0", "dap1", "das0", "das1", "pap0", "pap1")		##já que o arquivo já se chama ab, simplificamos os nomes das colunas

#Transformando as medidas dap e pap numa só: somatória dos valores dessas duas colunas por linha:
ab.final$dap1=rowSums(ab.final[,c("dap1","pap1")],na.rm=T)		##essa função rowSums aparece no help da sum(); essa indexação é muito importante de saber!
ab.final$dap0=rowSums(ab.final[,c("dap0","pap0")],na.rm=T)
head(ab.final)	##conferida a soma pelas plantas que tinham só pap

#Gerando um valor de dap final:
diam.final=sqrt(ab.final[,-1]*(4/pi))	##transformando a área basal em um valor de dap final, exceto pela coluna 1 de ab.final, que é a coluna de placa!
diam.final[diam.final==0]=NA 		## prestar atenção porque essa indexação é diferente(sem a vírgula!) é para matriz ou data.frame; aqui não queremos manter os zeros inseridos quando fizemos as somatórias, por isso, passamos tudo que é zero para NA

#Agora retornamos a coluna placa para esse novo objeto diam.final
diam.final$placa=ab.final$placa

#Só que a coluna placa foi adicionada no final e não precisamos mais das colunas de pap0 e pap1
diam.final=diam.final[,c(7,1:4)]	##reordenamos as colunas por indexação e ainda removemos as colunas 5 (pap0) e 6 (pap1)
diam.final		##ok, agora esse arquivo já tem os diâmetros finais para os dois tempos (dap e das)


### Preparando os dados para as transições ###

#Criando Status para o diam.final: precisamos adicionar colunas com os status0 e status1:

#Ao invés de buscar essas colunas no objeto teste, onde cada linha é um fuste, vamos criar as colunas status, pois são três apenas os estados possíveis.
	#Status possíveis: 
	#V=viva
	#M=morta
	#NE=não existe	status para as plantas recrutadas num tempo anterior ao recrutamento

#Primeiro criamos as duas novas colunas desses status com NA's:
diam.final$status0=NA
diam.final$status1=NA
diam.final	##conferindo...

#Definindo as condições para cada um dos status:

#Para o tempo 0:
diam.final$status0[is.na(diam.final$dap0) & is.na(diam.final$das0) ]="NE"	##na coluna status0 todos os indivíduos que tiverem dap0 e das0 iguais a "NA" atribuímos "NE" (ou seja, todos que não existiam no 0 mas existem no 1)
diam.final$status0[diam.final$dap0>0 | diam.final$das0>0]="V"	##na coluna status0 todos os indivíduos com dap0 ou das0 maior do que 0

#Para o tempo 1:
diam.final$status1[diam.final$status0=="V" & is.na(diam.final$dap1) & is.na(diam.final$das1)]="M"	##na coluna status1 todos os indivíduos que tem status0=V e que são NA no dap1 e no das1 atribuímos "M"
diam.final$status1[diam.final$dap1>0 | diam.final$das1>0]="V"	##por outro lado, todos em que o dap1 ou o das1 são maiores do que zero atribuímos "V"

#Uma questão importante: nunca fazer comandos redundantes, pois dessa forma não podemos enxergar possíveis erros. 
#Por isso para a atribuição de cada status0 e 1 na primeira expressão fazemos uma indexação e na segunda outra e não no segundo comando simplesmente mandamos chamar de outro status tudo que não estava contemplado no primeiro

#table(is.na(diam.final$dap1))	##apenas conferindo se havia mortas (quantos NA's no dap1 que não são das...)
#diam.final[is.na(diam.final$dap1),]	##quem eles são

diam.final$status0
diam.final$status1
diam.final

#placas.prob=diam.final[is.na(diam.final$status1),"placa"]	##apenas uns testes para os NA's que tinham aparecido estranhamente... mas já foi solucionado, só mantive aqui para guardar o comando match
#match(placas.prob,teste$placa)
#teste[match(placas.prob,teste$placa),]

diam.final[is.na(diam.final$status1),]	##checando os NA's após a correção mencionada acima, resta apenas a placa 4185, que está sem nenhum valor de medida

#Uma vez que os diâmetros finais estão prontos e os status também, podemos estabelecer as classes de tamanho:

### Estabelecendo as classes ### 

#Checagem dos mínimos e máximos:
min(diam.final$das0,na.rm=T)
max(diam.final$das0,na.rm=T)
hist(diam.final$das0)

#Criação de uma coluna para a classe de tamanho do tempo zero (primeiro com NA's)
diam.final$classe0=NA

#Atribuição de classes a essa coluna:
diam.final$classe0[diam.final$das0>0]="p"		##classe "plântula" para todos que tem das
diam.final$classe0[diam.final$dap0>=10 & diam.final$dap0<40]="j1"		##classe "juvenil 1"; notar que a indexação é feita em dois momentos!
diam.final$classe0[diam.final$dap0>=40 & diam.final$dap0<70]="j2"		##classe "juvenil 2"
diam.final$classe0[diam.final$dap0>=70 & diam.final$dap0<100]="a1"	##classe "adulto 1"
diam.final$classe0[diam.final$dap0>=100]="a2"		##classe "adulto 2"
diam.final

#Criação de uma coluna para a classe de tamanho do tempo 1 (primeiro com NA's)
diam.final$classe1=NA

diam.final$classe1[diam.final$das1>0]="p"
diam.final$classe1[diam.final$dap1>=10 & diam.final$dap1<40]="j1"
diam.final$classe1[diam.final$dap1>=40 & diam.final$dap1<70]="j2"
diam.final$classe1[diam.final$dap1>=70 & diam.final$dap1<100]="a1"
diam.final$classe1[diam.final$dap1>=100]="a2"
diam.final

### As transições, finalmente! ###

#Preciso contar quantos permaneceram na mesma classe e quantos mudaram:

#Começando pelos que permaneceram na mesma classe:

table(diam.final$classe0==diam.final$classe1)
#FALSE  TRUE 
#   9   181 	##só nove transições! 

#Mas são 192 observações, cadê as outras 2?
#Ah! são NA's nas classes, problemas a conferir na planilha original...

#Quais são então as transições que acontecem NESSE CONJUNTO DE DADOS?

diam.final[diam.final$classe0!=diam.final$classe1,]
# p -> j1
# j1 -> j2
# j2 -> a1
# a1 -> a2

table(diam.final$classe0)
table(diam.final$classe1)

#Achei que um comando desse tipo faria a contagem das combinações entre classes, mas não deu certo!
#xtabs(cbind(classe0, classe1) ~ ., data = diam.final))	##não funciona!

#Calculando as transições uma a uma:

# p -> p: permanência
diam.final[diam.final$classe0=="p" & diam.final$classe1=="p",]	##indexando linha e coluna ele mostra quais são as linhas de dados que satisfazem as condições
table(diam.final$classe0=="p" & diam.final$classe1=="p")	##fazendo um teste lógico combinado com o table ele tabula o número de verdadeiros e falsos
#FALSE  TRUE 
#  145    46
table(diam.final$classe0=="p")	##quantos indivíduos pertenciam à classe p no tempo 0
#FALSE  TRUE 
#  143    48 
tr.pp <- round(46/48,4)		##não consegui automatizar isso!
tr.pp
[1]  0.9583


####################TENTEI CRIAR UM COMANDO GENÉRICO, MAS NÃO DESCOBRI ######################################################################
#tr.pp <- table(table(diam.final$classe0=="p" & diam.final$classe1=="p")=="F")	##como faço para ele dizer o número de falsos?
#tr.pp <- round(table(table(diam.final$classe0=="p" & diam.final$classe1=="p")=="F")/192, 4)	##como faço para ele dizer o número de falsos?
#tr.pp
#dim(diam.final[1])	##como faço para dizer para ele pegar o número de linhas do objeto e utilizar no cálculo?
##############################################################################################################################################

# p -> j1: crescimento
diam.final[diam.final$classe0=="p" & diam.final$classe1=="j1",]
table(diam.final$classe0=="p" & diam.final$classe1=="j1")
#FALSE  TRUE 
#  189     2 
table(diam.final$classe0=="p")
#FALSE  TRUE 
#  143    48 
tr.pj1 <- round(2/48,4)
tr.pj1
#[1] 0.0417

# p -> j2: crescimento
diam.final[diam.final$classe0=="p" & diam.final$classe1=="j2",]
table(diam.final$classe0=="p" & diam.final$classe1=="j2")
# FALSE 
# 191 
tr.pj2 <- 0/48
tr.pj2
#[1]  0

# p -> a1: crescimento - é zero, mas preciso calcular aqui?
diam.final[diam.final$classe0=="p" & diam.final$classe1=="a1",]
table(diam.final$classe0=="p" & diam.final$classe1=="a1")
# FALSE 
# 191 
tr.pa1 <- 0/48
tr.pa1
#[1]  0

# p -> a2: crescimento - é zero, mas preciso calcular aqui?
diam.final[diam.final$classe0=="p" & diam.final$classe1=="a2",]
table(diam.final$classe0=="p" & diam.final$classe1=="a2")
# FALSE 
# 191 
tr.pa2 <- 0/48
tr.pa2
#[1]  0

# j1 -> p: retração
diam.final[diam.final$classe0=="j1" & diam.final$classe1=="p",]
table(diam.final$classe0=="j1" & diam.final$classe1=="p")
# FALSE 
# 191 
table(diam.final$classe0=="j1")
#FALSE  TRUE 
#  125    66 
tr.j1p <- 0/66
tr.j1p
#[1]  0

# j1 -> j1: permanência
diam.final[diam.final$classe0=="j1" & diam.final$classe1=="j1",]
table(diam.final$classe0=="j1" & diam.final$classe1=="j1")
#FALSE  TRUE 
# 126    65 
tr.j1j1 <- round(65/66,4)
tr.j1j1
#[1] 0.9848


# j1 -> j2: crescimento
diam.final[diam.final$classe0=="j1" & diam.final$classe1=="j2",]
table(diam.final$classe0=="j1" & diam.final$classe1=="j2")
#FALSE  TRUE 
#  190     1 
tr.j1j2 <- round(1/66,4)
tr.j1j2
#[1] 0.0152


# j1 -> a1: crescimento
diam.final[diam.final$classe0=="j1" & diam.final$classe1=="a1",]
table(diam.final$classe0=="j1" & diam.final$classe1=="a1")
#FALSE 
#  191 
tr.j1a1 <- round(0/66,4)
tr.j1a1
#[1] 0


# j1 -> a2: crescimento
diam.final[diam.final$classe0=="j1" & diam.final$classe1=="a2",]
table(diam.final$classe0=="j1" & diam.final$classe1=="a2")
#FALSE 
#  191 
tr.j1a2 <- round(0/66,4)
tr.j1a2
#[1] 0


# j2 -> p: retração
diam.final[diam.final$classe0=="j2" & diam.final$classe1=="p",]
table(diam.final$classe0=="j2" & diam.final$classe1=="p")
#FALSE 
#  191 
table(diam.final$classe0=="j2")
#FALSE  TRUE 
#  161    30 
tr.j2p <- round(0/30,4)
tr.j2p
#[1] 0


# j2 -> j1: retração
diam.final[diam.final$classe0=="j2" & diam.final$classe1=="j1",]
table(diam.final$classe0=="j2" & diam.final$classe1=="j1")
#FALSE 
#  191 
tr.j2j1 <- round(0/30,4)
tr.j2j1
#[1] 0


# j2 -> j2: permanência
diam.final[diam.final$classe0=="j2" & diam.final$classe1=="j2",]
table(diam.final$classe0=="j2" & diam.final$classe1=="j2")
#FALSE  TRUE 
#  163    28 
tr.j2j2 <- round(28/30,4)
tr.j2j2
#[1]  0.9333


# j2 -> a1: crescimento
diam.final[diam.final$classe0=="j2" & diam.final$classe1=="a1",]
table(diam.final$classe0=="j2" & diam.final$classe1=="a1")
#FALSE  TRUE 
#  189     2
tr.j2a1 <- round(2/30,4)
tr.j2a1
#[1]   0.0667


# j2 -> a2: crescimento
diam.final[diam.final$classe0=="j2" & diam.final$classe1=="a2",]
table(diam.final$classe0=="j2" & diam.final$classe1=="a2")
#FALSE  
#  191  
tr.j2a2 <- round(0/30,4)
tr.j2a2
#[1] 0


# a1 -> p: retração
diam.final[diam.final$classe0=="a1" & diam.final$classe1=="p",]
table(diam.final$classe0=="a1" & diam.final$classe1=="p")
#FALSE  
#  190
table(diam.final$classe0=="a1")
#FALSE  TRUE 
#  170    21 
tr.a1p <- round(0/21,4)
tr.a1p
#[1] 0


# a1 -> j1: retração
diam.final[diam.final$classe0=="a1" & diam.final$classe1=="j1",]
table(diam.final$classe0=="a1" & diam.final$classe1=="j1")
#FALSE  
#  190
tr.a1j1 <- round(0/21,4)
tr.a1j1
#[1] 0


# a1 -> j2: retração
diam.final[diam.final$classe0=="a1" & diam.final$classe1=="j2",]
table(diam.final$classe0=="a1" & diam.final$classe1=="j2")
#FALSE  
#  190
tr.a1j2 <- round(0/21,4)
tr.a1j2
#[1] 0


# a1 -> a1: permanência
diam.final[diam.final$classe0=="a1" & diam.final$classe1=="a1",]
table(diam.final$classe0=="a1" & diam.final$classe1=="a1")
#FALSE  TRUE 
#  174    16 
tr.a1a1 <- round(16/21,4)
tr.a1a1
#[1] 0.7619


# a1 -> a2: crescimento
diam.final[diam.final$classe0=="a1" & diam.final$classe1=="a2",]
table(diam.final$classe0=="a1" & diam.final$classe1=="a2")
#FALSE  TRUE 
#  186    4 
tr.a1a2 <- round(4/21,4)
tr.a1a2
#[1]  0.1905

##AQUI NÃO ESTÁ SOMANDO 1!ACHO QUE SÃO OS NA'S

# a2 -> p: retração
diam.final[diam.final$classe0=="a2" & diam.final$classe1=="p",]
table(diam.final$classe0=="a2" & diam.final$classe1=="p")
#FALSE  
#  191
table(diam.final$classe0=="a2")
#FALSE  TRUE 
#  165    26 
tr.a2p <- round(0/26,4)
tr.a2p
#[1] 0


# a2 -> j1: retração
diam.final[diam.final$classe0=="a2" & diam.final$classe1=="j1",]
table(diam.final$classe0=="a2" & diam.final$classe1=="j1")
#FALSE  
#  191
tr.a2j1 <- round(0/26,4)
tr.a2j1
#[1] 0


# a2 -> j2: retração
diam.final[diam.final$classe0=="a2" & diam.final$classe1=="j2",]
table(diam.final$classe0=="a2" & diam.final$classe1=="j2")
#FALSE  
#  191
tr.a2j2 <- round(0/26,4)
tr.a2j2
#[1] 0


# a2 -> a1: retração
diam.final[diam.final$classe0=="a2" & diam.final$classe1=="a1",]
table(diam.final$classe0=="a2" & diam.final$classe1=="a1")
#FALSE  
#  191
tr.a2a1 <- round(0/26,4)
tr.a2a1
#[1] 0


# a2 -> a2: permanência
diam.final[diam.final$classe0=="a2" & diam.final$classe1=="a2",]
table(diam.final$classe0=="a2" & diam.final$classe1=="a2")
#FALSE  TRUE 
#  165    26 
tr.a2a2 <- round(26/26,4)
tr.a2a2
#[1] 1

### AGORA JUNTANDO TUDO NUMA MATRIZ ###

dados <- c(tr.pp, tr.pj1, tr.pj2, tr.pa1, tr.pa2, tr.j1p, tr.j1j1, tr.j1j2, tr.j1a1, tr.j1a2, tr.j2p, tr.j2j1, tr.j2j2, tr.j2a1, tr.j2a2, tr.a1p, tr.a1j1, tr.a1j2, tr.a1a1, tr.a1a2, tr.a2p, tr.a2j1, tr.a2j2, tr.a2a1, tr.a2a2)
trans.matrix <- matrix(data = dados, nrow = 5, ncol = 5, byrow = FALSE)
trans.matrix
colnames (trans.matrix) <- c("p", "j1", "j2", "a1", "a2")
rownames (trans.matrix) <- c("p", "j1", "j2", "a1", "a2")
trans.matrix


</code>