======Seleção de modelos - Alometria======

===== Agosto 09 =====

==== Seleção e decisões finais ====
  * A princípio decidimos verificar se havia efeito do grupo (leg e não-leg). Essa diferença existe e esta demonstrada no script abaixo. Mas, depois mudamos de idéia (ver resultados e explicações a seguir)

{{script_models_alo.txt|Script seleção de modelos GRUPO}} 

  * No final, decidimos agrupar leg e não leg (apesar do efeito de grupo encontrado) porque para comparar temos que ter uma área igual no entorno de leg e não-leg!!!!
  * Os três modelos contrastados foram: 
    - Raio copa não varia em função do DAP
    - Raio copa varia linearmente em função do DAP
    - Raio copa varia exponencialmente em função do DAP
  * Analisando os valores de AIC, o modelo mais plausível foi o 3, e a equação encontrada nele foi: **0.047*(DAP)^0.8**.
  * Script dessa seleção:
 
{{script_selecao_alometria_abril2010.txt| Seleção de modelos DAP x COPA}}

==== Coleta dos dados ====
  * Amostra estratificada e balanceada
  * Dados de DAP e altura
  * Quatro medidas de raio copa, em cruz a partir do tronco, tirei média depois.

{{  graf_sep.jpeg?800  |}} 

{{copa_alo_uso.csv|Dados brutos finais}}


=====Piloto Alometria Julho 09=====

====Coleta====

  * Coletei no campo dados de diâmetro maior e diâmetro menor da copa (cm)
  * Calculei a área considerando uma elipse -> pi * (DMaior/2) * (DMenor/2)
  * Dúvida: chego no raio considerando que essa área é de uma circunferência?
  * Acho que to errando nas transformações, pois gero raios de copa gigantes...
{{dados_maio09.csv|Dados brutos piloto}}



==== Seleção modelos =====
Escolhi as árvores estratificadas, entretanto, em campo algumas não leg estavam mortas, outras emaranhadas com outra copa impossibilitando a visualização e algumas não foram encontradas. Algumas árvores estavam em parcelas que ainda não foram recenseadas e o acesso estava muito difícil.
Acima de 450 de dap, a grande maioria é balizia e do outro grupo há alguns Calophyllum e Ocotea.

{{gráf.copaxdap_grupo.jpeg?500|}}

Testando efeito do grupo
  - Modelo 1 = copa ~ dap
  - Modelo 2 = copa ~ dap + grupo
  - Modelo 3 = copa ~ dap * grupo
Comparando com anova e AIC, há interação... 


  - Modelo linear para legs
    * Regressão simples: **Raio.copa.leg = 0.49 + 0.01*(Dap.leg)**... significativa, com r² = 0.9.  
  - Modelos para não legs
    * Regressão simples: **Raio.copa.nleg = 1.3 + 0.004* Dap.nleg**... significativa com r² = 0.5
    * M-M: Raio.copa.nleg = 2.92*Dap.nleg/(56.6 + Dap.nleg)
Os AICs foram similares, fiquei com a regressão simples. Além disso, o erro da parâmetro b na MM foi enorme (26!).

==== Arquivos R====

{{alo.rdata|Arquivo RData com modelos de alometria julho 09}}    

{{script_alo_juntas.r|Script 1ª análise(leg + nleg) julho09}}

{{script_alo_legxnleg.r|Script 2ª análise (efeito do grupo) julho09}}

==== Orientações gráfico ====

  plot(rm.copa~dap, col=c("red", "blue")[as.factor(copa$Grupo)],pch=19, xlab= "DAP (cm)", ylab= "raio da copa (m)")
  ## rm.copa é o raio em metros da copa e o dap é o dap em cm. Copa$grupo é a coluna de dados com o  fator leg x nleg 
  
   locator(1)->locxy
  ## aqui vc clica  no ponto xy no gráfico onde deseja colocar a legenda e guarda no objeto locxy.
  
  legend(locxy, legend=c("leguminosae", "outras familias"), col=c("red", "blue"), pch=c(19,19))
  ## inclui a legenda na posição locxy

{{graf_dapxcopa1.jpg?500|}}


==== Legs juntas ====
Como eram poucos dados plotei todos leg e não leg juntas...  \

{{gráfico.jpeg?500|}}  \

Comparei 5 modelos:  \
1)O mais simples : 
<code>
copa.vero1<-function(m=2.542,s=1.96)
	{
     -sum(stats::dnorm(alo$Raio.copa.m, m, s, log=TRUE))
	}
copa1<-mle(copa.vero1)
</code>
2) Com a função original: 
<code>
copa.vero2<-function(a=0,b=1,s=1.96)
	{
	m = a+b*((0.1*alo$Dap.mm)^0.5/2)
	-sum(stats::dnorm(x=alo$Raio.copa.m,mean=m,sd=s,log=T))
	}
</code>
copa2<-mle(copa.vero2)
3) Com outra função, sem elevar a 0.5: 
<code>
copa.vero3<-function(a=0,b=1,s=1.96)
	{
	m = a+b*((0.1*alo$Dap.mm)/2)
	-sum(stats::dnorm(x=alo$Raio.copa.m,mean=m,sd=s,log=T))
	}
</code>
copa3<-mle(copa.vero3)
4) Elevando ao quadrado sem dividir: 
<code>
copa.vero4<-function(a=0,b=1,s=1.96)
	{
	m = a+b*(0.1*alo$Dap.mm)^2
	-sum(stats::dnorm(x=alo$Raio.copa.m,mean=m,sd=s,log=T))
	}
</code>
copa4<-mle(copa.vero4)

5) Com média e desvio variando pela mesma função: 
<code>
copa.vero5<-function(a=1.43,b=0.002,s1=1,s2=2)
	{
	m = a+b*((0.1*alo$Dap.mm)/2)
	s = s1+s2*((0.1*alo$Dap.mm)/2)
	-sum(stats::dnorm(x=alo$Raio.copa.m,mean=m,sd=s,log=T))
	}
</code>
copa5<-mle(copa.vero5)  

| |AIC|
|copa1|111.70791|
|copa2|71.15744|
|copa3|68.34604|
|copa4|77.40857|
|copa5|63.01979|

* Como a curiosidade é grande rodei alguns nulos alterando a função r.copa para o modelo 3 e o resultado foi:
      * Para Abarbr, 259 médias ficaram inferiores ao valor observado. Com a função anterior esse valor era 0.
      * Balipe continuo com 0, como possui bastante dbh grande acredito que a função possa estar errando aí...
      * Ormoar teve 8 valores inferiores, com a função anterior era 1 (sem muita diferença então)
====Comentário Ale=====
Continuo não entendendo porque os modelos tem parâmetros como "0.01*" e "/2". Não fazem o minimo sentido. Esses são os parâmetros que vc. quer estimar... eles estão incluídos no "a" e "b" da regressão linear...

Acho que devemos nos limitar nesse momento apenas ao lm.
Minhas sugestões talvez não tenham ficado claras. 

1. passo:
Confrontar o modelo linear com o fator "leg" e "nleg", com o modelo linear se esse fator, se não me engano vc. chamou "Grupo".
  * mod1=lm(copa~dap*grupo)
  * mod2=lm(copa~dap)
use: summary(); anova() para ver os modelos separadamente e depois anova ou AICtab para confrontá-los. Se no modelo um há interação significa que leg e nleg devem ser tratados diferentemente na hora de estimar a copa. Caso o modelo só com o dap seja melhor, então podemos fazer um modelo único.

2.passo: decidido qual o modelo base, devemos estimar os parâmetros (leg/nleg junto ou isolado dependendo do passo anterior)
  * modelo linear simples e com desvio padrão variando 

   copa~a+ b* dap

3. estimado a e b voltamos ao modelo nulo para ajustá-lo

==== Comentários Ale ====

  * mesmo para um piloto os dados estão mal coletados. Não há a distribuição de diâmetros estratificada como combinamos e tb. não há correspondência leg x nleg, como vc. afirmou! Para evitar essas surpresas deve olhar os dados ainda quando está em campo para fechar lacunas de coleta. Com esses dados não há como chegar a um modelo legal de predição, mas ao menos podemos fazer um ensaio para esse relatório. Assim os códigos ficam prontos para serem usados quando vc. coletar mais dados.
 
  * a partir desse relatórios não estamos mais na fase de pilotos e tentivas. Agora é vida real e vc. precisa se conscientizar que os prazos estão contra vc.... infelizmente não há mais tempo para ficarmos contemplando e pirando em modelos. Também não há mais tempo para vc. ficar quebrando a cabeça com códigos. QQ travada que der deve me procurar imediatamente para resolvermos. 

  * até agora há mais esforço meu que seu no trabalho intelectual de definição dos modelos nulos. Eu gosto da brincadeira, mas o brinquedo não é meu!  A partir de agora vc. assume o brinquedo e eu apenas facilito sua vida com os códigos e ajudo apontando caminhos. A decisão de até onde sua monografia irá chegar é sua.

  * nesse sentido estabeleço o que julto ser o mínimo para ser atingido (ao final da tese):
    - alometria da copa x diâmetro para melhorar o modelo nulo (n min: 100 arvores legs e 100 nleg), estratificado e balanceado 
    - testar o modelo nulo de índices para todas as legs (usando os dados do recenso!)
    - testar para jovens com DAP<1cm leg x nleg (levantar os jovens em campo), n mínimo de parcelas: 30 abaixo de legs e 30 abaixo de nlegs

  * até o dia da entrega do seu relatório vc. deve atualizar o que fez diariamente no wiki. O que vamos fazer é:
     - seleção do modelo de predição de raio de copa como os dados que tem.
     - ajustar o modelo nulo com essa alometria e rodar para todas as legs.
     - as análises devem ser feitas essa semana para vc. poder escrever o relatório na semana que vem.
  
  * depois do seu relatório vc. deve atualizar mensalmente o wiki e prever um boa e longa estada em campo para:
     - levantar os dados de alometria copa x dap
     - levantar dados de jovens dap< 1cm  abaixo de leg e nleg
     - agilizar a entrada de dados da parcela permanente para estar apta a usar os dados ao final do seu trabalho
  * durante esse tempo deve focar em leitura específica sobre o seu trabalho e mais ainda sobre no artigo que irá apresentar na sua qualificação. Quando terminar o relatório decidimos qual será o artigo...