Métodos de Calibração de
Modelos hidrológicos
Carlos Ruberto Fragoso Júnior
Marllus Gustavo Ferreira Passos das Neves
Sumário





Conceito básicos
 O que é calibração?
 Problemas comuns na calibração de modelos hidrológicos
 Ciclo da calibração
Métodos de calibração
Função objetivo
Técnicas numéricas
 Busca aleatória
 Técnicas iterativas;
 Busca direta;
 Técnicas de otimização global;
 Algoritmos genéticos
Critérios de parada
O que é calibração?

Procura de valores dos parâmetros de um modelo matemático
que resultem em uma boa concordância entre dados
observados e calculados
O erro é minimizado!
Calibração - Otimização

Otimização  encontrar o mínimo ou o
máximo de uma função
160
140
120
100
80
60
40
20
0
0
10
20
30
Calibração - Otimização

Otimização  calibração: quando a função a
ser minimizada representa a discrepância entre
os dados calculados por um modelo e os dados
medidos

Calibração  o modelo “quer ser” o sistema
que ele vai simular
Problemas comuns em modelos hidrológicos

Encontrar um conjunto ótimo de parâmetros que
ajusta um evento de cheia ou uma série de
vazões ou vazões mínimas

Exemplo  encontrar o coeficiente do
reservatório linear simples que ajusta
adequadamente uma recessão de um
hidrograma
Q(t+Dt) = Q(t) . exp(-Dt/k)
Exemplo
Exemplo
3:50
Exemplo
Q(t+Dt) = Q(t) . exp(-Dt/k)
Primeiro teste: k = 20
Problemas na calibração de modelos hidrológicos




Modelos hidrológicos geralmente tem muitos
parâmetros
Não lineares
Técnicas de otimização automáticas
Usar Funções Objetivo
Ciclo da calibração
Critérios para mudança
dos parâmetros
Critérios para um “bom
ajuste” (Função objetivo)
Rodar o modelo
Novos valores de
parâmetros
Critérios de parada
Verificar o erro
ou FO
Métodos de calibração
Métodos de
calibração
Tentativa e erro
(Manual)
Técnicas numéricas
Ajusta os parâmetros
manualmente
baseado nos
resultados
Usa algoritmos
numéricos para
encontrar um conjunto
de parâmetros
ótimo
Aleatório
Assume faixa de
probabilidade para
cada parâmetro
Funções Objetivo (FO)

Medida do erro – objetivo é minimizar a FO
Diferentes funções objetivo







Somatório dos erros: compensação de erros
Somatório do módulo dos erros
Somatório dos erros ao quadrado
Somatório de erros relativos
Somatório dos desvios dos inversos da vazão
Erro de volume (bias)
Coeficiente de eficiência de Nash-Sutcliffe
70
Observada
60
simulada
50
vazões

40
30
20
10
0
0
20
40
60
80
tempo
100
120
140
Funções Objetivo
N
F1   (QOi  QCi )
2
i 1
N
F2   QO i  QC i
i 1
N
1
1 2
F3   (

)
QO
QC
i
i
i 1
N
F4   (
i 1
N
QOi  QCi 2
)
QOi
Função quadrática  dá um peso maior para
as vazões maiores  a tendência será
ajustar melhor as máximas (série contínua)
Função módulo  também dá maior peso às
vazões maiores e é pouco utilizada
(superfície não comportada)
Função para mínimos  dá um peso maior
para as vazões mínimas  inverte a
prioridade de F1
Função relativa  procura retirar o peso
excessivo dado aos valores absolutos de F1
2
F5   (logQOi  logQCi )
i1
Exemplo
N
F1   (QOi  QCi )
2
i 1
N
F2   QO i  QC i
i 1
N
F3   (
i 1
3:50
1
1 2

)
QOi QCi
Funções Objetivo

Raiz do Erro Médio Quadrado (RSME)
RMSE 


n
i 1
( X obs,i  X mo del ,i ) 2
n
Raiz do Erro Médio Quadrado Normalizado (NRSME)
RMSE
NRMSE 
X obs,max  X obs,min
Quanto mais perto de zero
melhor
RMSE
NRMSE 
X obs
Funções Objetivo

Coeficiente de correlação de Pearson
r



n
i 1
n
i 1
( xi  x)  ( yi  y )
( xi  x) 
2

n
i 1
( yi  y ) 2
Coeficiente de eficiência de Nash-Sutcliffe

E  1

n
i 1
n
( X obs ,i  X mo del ) 2
i 1
( X obs ,i  X obs ) 2
Quanto mais perto de 1 melhor
Funções Objetivo - comparação

Caso
Ideal
E=1
RMSE = 0
Funções Objetivo - comparação

Caso 1
Máximas
próximas,
mínimas
dispersas, tenta
acompanhar o
formato
3,20
3,00
2,80
2,60
2,40
2,20
2,00
1,80
1,60
1,40
1,20
1,00
0,80
0,60
0,40
0,20
0,00
ymodel
yobs
Média obs
3,50
R² = 0,8668
0 0,53,00
1 1,5 2 2,5 3 3,5 4 4,5 5 5,5 6 6,5 7 7,5 8 8,5 9 9,5 1010,51111,5
2,50
E = 0,791
RMSE = 0,318
2,00
1,50
1,00
0,50
0,00
0,00
0,50
1,00
1,50
2,00
2,50
3,00
3,50
Funções Objetivo - comparação

Caso 2
Mínimas
próximas,
máximas
dispersas, tenta
acompanhar o
formato
3,20
3,00
2,80
2,60
2,40
2,20
2,00
1,80
1,60
1,40
1,20
1,00
0,80
0,60
0,40
0,20
0,00
ymodel
yobs
Média obs
3,50
0 0,5 1 1,5 2 2,5 3 3,5 4 4,5 5 5,5 6 6,5 7 7,5 8 8,5 9 9,5 1010,51111,5
3,00
R² = 0,8633
2,50
E = 0,843
RMSE = 0,276
2,00
1,50
1,00
0,50
0,00
0,00
0,50
1,00
1,50
2,00
2,50
3,00
3,50
Funções Objetivo - comparação

Caso 3
Máximas
próximas,
mínimas
dispersas, não
acompanha o
formato.
Observar os
pontos em torno
da média
E = 0,263
RMSE = 0,598
3,20
3,00
2,80
2,60
2,40
2,20
2,00
1,80
1,60
1,40
1,20
1,00
0,80
0,60
0,40
0,20
0,00
ymodel
yobs
Média obs
3,50
3,00
R² = 0,4031
0 0,5 1 1,5 2 2,5 3 3,5 4 4,5 5 5,5 6 6,5 7 7,5 8 8,5 9 9,5 1010,51111,5
2,50
2,00
1,50
1,00
0,50
0,00
0,00
0,50
1,00
1,50
2,00
2,50
3,00
3,50
Funções Objetivo - comparação

Caso 4
Mínimas
próximas,
máximas
dispersas, não
acompanha o
formato.
Observar os
pontos em torno
da média
E = 0,541
RMSE = 0,471
3,20
3,00
2,80
2,60
2,40
2,20
2,00
1,80
1,60
1,40
1,20
1,00
0,80
0,60
0,40
0,20
0,00
ymodel
yobs
Média obs
3,50
0 0,5 1 1,5 2 2,5 3 3,5 4 4,5 5 5,5 6 6,5 7 7,5 8 8,5 9 9,5 1010,51111,5
3,00
2,50
2,00
1,50
1,00
R² = 0,6427
0,50
0,00
0,00
0,50
1,00
1,50
2,00
2,50
3,00
3,50
Funções Objetivo - comparação

Caso 5
Longe dos
máximos e dos
mínimos
3,20
3,00
2,80
2,60
2,40
2,20
2,00
1,80
1,60
1,40
1,20
1,00
0,80
0,60
0,40
0,20
0,00
ymodel
yobs
Média obs
3,50
3,00
0 0,5 1 1,5 2 2,5 3 3,5 4 4,5 5 5,5 6 6,5 7 7,5 8 8,5 9 9,5 1010,51111,5
2,50
E = 0,573
RMSE = 0,455
2,00
1,50
1,00
R² = 0,5997
0,50
0,00
0,00
0,50
1,00
1,50
2,00
2,50
3,00
3,50
Técnicas de otimização


Cálculo analítico
Técnicas numéricas



Busca aleatória
Busca direta
Algoritmos genéticos
Cálculo analítico
Encontrar pontos da função em que a derivada é
zero  pode ser rápido, mais elegante, mas depende de
alguns requisitos via de regra incomuns nas F.O. de
calibração de um modelo chuva-vazão (funções de picos
múltiplos, funções descontínuas, ausência da forma
analítica da função)
160
140
F(x)  a  b.x  c.x 2
120
100
80
60
40
dF
0
dx
20
0
0
10
20
30
Cálculo analítico
2 variáveis: F.O. analítica F: D  R, D C R2 fechado e
limitado, f diferenciável e de classe C2  f assume um valor
máximo e um mínimo em D (teorema de Weierstrass)

Condição necessária  existam pontos críticos  F  0
 2F
Condição suficiente  incluir o estudo
x2
da função Hessiana num ponto interior H(x, y)   2F
(xo,yo)
yx
(a)
(b)
(c)
(d)
Se H(xo,yo)
Se H(xo,yo)
Se H(xo,yo)
Se H(xo,yo)
2F
xy
2F
y2
> 0 e Fxx(xo,yo) > 0  (xo,yo) é ponto de mínimo local
> 0 e Fxx(xo,yo) < 0  (xo,yo) é ponto de máximo local
< 0  ponto de sela
= 0  nada se afirma
Exemplo
Determine o mínimo da função
F
 3x2  3  0
x
F
 3y2  3  0
y
F(x, y)  x3  y3  3x  3y  4
Pontos críticos: (1,1), (1,-1), (-1,1), (-1,-1)
Os pontos (1,1) e (-1,-1) são candidatos
a mínimo e máximo, respectivamente:
pelo cálculo das derivadas segundas
2F
2F
(-1,-1)  - 6  2 (-1,-1)
2
x
y
2F
2F
(1,1)  6  2 (1,1)
2
x
y
6x 0
H(x, y) 
 36xy As candidaturas se confirmam
0 6y
Características das Técnicas Numéricas
Técnicas automáticas
(Nash e Sutcliffe, 1970):
Envolve sucessivas variações dos valores
dos parâmetros de acordo com alguma
regra ou padrão pré-concebidos de
incrementos que levam em conta os
resultados dos passos anteriores e em
particular se a variação melhorou ou não o
ajuste.
Características das Técnicas Numéricas
Definição do ponto de partida: critério para
inicializar o processo de tentativa  em geral
depende mais do problema em questão do que do
método
Direção de pesquisa: identifica o vetor no qual serão
realizadas as alterações das variáveis
Espaçamento de cada tentativa: variação que
ocorrerá na direção de pesquisa a cada tentativa
Critérios de parada: definição dos critérios para
aceitar uma determinada solução como o ótimo de
uma função
Técnicas numéricas - Busca Aleatória
Conceito de probabilidade
Mais simples  Geram-se N valores dos parâmetros,
calculam-se as F.O. e escolhe-se o menor valor de F.O
Vantagens: funções
descontínuas; picos
múltiplos
Desvantagens:
demorado; não existe
garantia de atingir o
ponto ótimo global
160
“Ótimo”
140
120
100
80
60
40
20
0
0
10
20
30
Técnicas numéricas - Busca direta
5
4.5
4
3.5
3
2.5
2
1.5
1
0.5
0
0

0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
Estratégia de caminhar “morro acima”
Função objetivo: F(x1,x2)
5
4.5
4
3.5
x2
3
2.5
2
1.5
Máximo local
1
Máximo global
0.5
0
0
0.5
1
1.5
2
2.5
x1
3
3.5
4
4.5
5
Início: ponto coordenadas (parâmetros) aleatórias
5
4.5
X2 = valor
aleatório entre
0e5
4
3.5
3
2.5
2
1.5
1
0.5
0
0
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
X1 = valor aleatório entre 0 e 5
Determina direção de busca: exemplo x2 = x2 + 0,3; x1 = x1
5
4.5
4
3.5
3
2.5
2
1.5
1
0.5
0
0
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
Função objetivo melhorou? Não, então tenta no outro sentido
5
4.5
4
3.5
3
2.5
2
1.5
1
0.5
0
0
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
F.O. melhorou? Sim, então continua no mesmo sentido
5
4.5
4
3.5
3
2.5
2
1.5
1
0.5
0
0
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
F.O. melhorou? Sim, então continua no mesmo sentido
5
4.5
4
3.5
3
2.5
2
1.5
1
0.5
0
0
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
F.O. melhorou? Sim, então continua no mesmo sentido
5
4.5
4
3.5
3
2.5
2
1.5
1
0.5
0
0
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
F.O. melhorou? Não, então volta para o ponto anterior...
...e muda a direção de busca.
5
4.5
4
3.5
3
2.5
2
1.5
1
0.5
0
0
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
F.O. melhorou? Sim, então continua no mesmo sentido
E assim segue até encontrar um ponto em que não existe
direção de busca que melhore o valor da F.O.
5
4.5
4
3.5
3
2.5
2
1.5
1
0.5
0
0
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
Método unidirecional ou univariacional
É o método anterior
1. Direção de pesquisa paralela aos
eixos;
2. Pesquisa em cada direção:
espaçamento constante ou variável
3. Critério de parada
desvantagens: (ao lado)
Método da rotação das coordenadas (Rosenbrook)
1º ciclo igual ao univariacional
 2º ciclo com rotação
 2 alternativas para pesquisa em cada direção: método
original que alterna a pesquisa de cada direção em cada
tentativa

Primeiro ciclo direção x1
x12  x10  S1  14; y2 = 88
x13  x11  S1  14  1,5x1  12,5; y3 = 63,25
x14  12,5  1,5x1,5  10,25 ; y4 = 34,56
x15
x16
Primeiro ciclo direção x2
x12  4  1  3; y8 = -12,01
x 22  3  1,5.1,0  1,5; y9 = -33,02
 10,25  1,5x 2,25  6,88; y5 = 10,52
x 32  1,5  1,5.1,5  0,75; y10 = -47,67
 6,88  1,5x3,375  1,82; y6 = 17,16
x 42  0,75  1,5.2,25  4,13; y11= -31,67
Rosenbrock: Método um pouco mais eficiente
5
4.5
4
3.5
3
2.5
2
1.5
1
0.5
0
0
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
Direção de busca é a que potencialmente dará maior incremento da FO
Limitação da busca direta: Ótimos locais
5
4.5
4
Região que
atrai solução
para o ótimo
local
3.5
3
2.5
2
1.5
1
0.5
0
0
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
Tentativa de contornar problema: Busca direta com
inicialização múltipla
5
4.5
Várias
tentativas;
espera se que
o ótimo global
seja a melhor
solução
testada.
4
3.5
3
2.5
2
1.5
1
0.5
0
0
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
Problema: Ineficiente e ineficaz quando a FO tem muitos ótimos locais
Técnicas numéricas – Busca direta

Busca direta (Rosenbrock e cia.)


vantagens: funções descontínuas; otimização
por simulação (funções que não podem ser
expressas analiticamente - calibração de
modelos)
desvantagens: funções com picos múltiplos
Técnicas numéricas – Busca direta
Técnicas numéricas – Algoritmos genéticos (AGs)
Pais mais adaptados têm maior probabilidade de gerar filhos
(sobrevivência do mais apto = seleção natural)
Darwin
Técnicas numéricas – Algoritmos genéticos (AGs)
• Conceitos da biologia (população, reprodução e gerações)
• Filhos são semelhantes aos pais (não completamente iguais)
• Pais mais “adaptados”  maior probabilidade de reproduzir


Natureza  indivíduos mais adaptados (com mais
“aptidão”)  maior probabilidade de sobreviver até chegar
à fase reprodutiva e participar do processo de
reprodução.
Algoritmo: pontos com maior F.O. (aptidão) têm maior
probabilidade de serem escolhidos para participar dos
complexos
Técnicas numéricas – Algoritmos genéticos (AGs)
Início
Inicialização da população
Cálculo da aptidão
Solução
encontrada?
Nova
população
Seleção
Reprodução
Mutação
Fim
Algoritmos genéticos (AGs)
Gera-se uma população
inicial (possíveis soluções)
População avaliada 
recebe uma nota ou
aptidão (qualidade da
solução)
Mais aptos selecionados e
os menos aptos
descartados
Os selecionados podem
sofrer operações de
crossover e mutação 
geram filhos  1ª geração
Nova avaliação  processo
repetido até se tornar
satisfatório
Algoritmo genético “puro”
1 - gera população (pontos aleatórios)
160
140
120
100
80
60
40
20
0
0
10
20
30
Algoritmo genético “puro”
2 - escolhe pontos para participar do processo de
“reprodução” (aqueles com melhor aptidão tem maior
probabilidade de escolha
160
140
120
100
80
60
40
20
0
0
10
20
30
Aptidão  pode ser igualada à F.O., mas via de regra, é
função da F.O. (evitar valores negativos ou convergência
prematura  extremos locais)
Algoritmo genético “puro”
3 - Exemplo de reprodução: escolhidos dois pontos  PAIS
160
140
120
100
80
60
40
20
0
0
10
Xa=8
20
Xb=19
30
 Binário:
Xa=01000
Xb=10011
Selecionados com probabilidade “proporcional à aptidão” 
mais aptidão, mais chance de serem escolhidos
Algoritmo genético “puro”
Genética: filhos “recebem” cromossomos dos
pais
01000
10011
01011
10000
Filhos:
É determinado um (ou mais) ponto de “corte”
(aleatório)
Filho 1: parte dos “cromossomos” do pai e parte da mãe
Filho 2: outra parte dos “cromossomos” do pai e parte da mãe
Xa=01011 = 11
Xb=10000 = 16
Algoritmo genético “puro”
Na operação de crossover, os pais transmitem aos filhos
algumas características
160
140
120
pais
100
80
60
40
filhos
20
0
0
10
20
30
Na mutação (a seguir), surgem novas características próprias
dos filhos  transmiti-las para outra geração
Algoritmo genético “puro”
Mutação: evento de baixa probabilidade
Genética: filhos “recebem”
cromossomos dos pais
01000
10011
É determinado um (ou mais) ponto de “corte”
(aleatório)
01011
10100
Filho 1: sem mutação
Filhos:
Filho 2: mutação
Xa=01011 = 11
Xb=10100 = 20
Algoritmo genético “puro”
Reprodução de todos os pontos escolhidos resulta na nova
geração
160
140
120
100
80
60
40
20
0
0
10
20
30
Mutação  melhora a diversidade na população, mas pode não
preservar a informação original  deve haver uma taxa
pequena (0,1 a 5%), mas suficiente
Algoritmo genético “puro”
Depois de algumas gerações
160
140
120
100
80
60
40
20
0
0
10
20
30
Algoritmo genético “puro”
Algumas desvantagens do algoritmo genético puro
• Números binários
• Transformação de variáveis de base decimal para binária
-0,05
+180,3
Variável Y
0000000000
1111111111
Usando 10 bits; Resolução = 0,176
decimal
Algoritmo genético “puro”
Algumas vantagens do algoritmo genético puro
Otimização com números inteiros
Diâmetros comerciais
Algoritmo genético “puro”


Representação binária
 Historicamente importante  tradicional
 Fácil de utilizar e manipular e simples de
analisar teoricamente
 Problema  com parâmetros contínuos e o
usuário desejando trabalhar com boa
precisão numérica  muito armazenamento
 convergência lenta do algoritmo
Representação real  melhor compreendida
 Há diversos operadores para n° real
(crossover e mutação)
Evolução de complexos misturados
(Shuffled complex evolution)


SCE - UA
Usa técnicas de





busca aleatória
algoritmos genéticos
simplex (Nelder e Mead)
Proposto por Duan, Gupta e
Sorooshian (U. Arizona)
Descrito no livro Sistemas
Inteligentes da ABRH
SCE-UA
Espaço de busca
Cria-se uma população de s pontos
de forma aleatória
População dividida em
algumas comunidades
(complexos ou conjuntos)
Evoluem de forma
independente  o espaço
de busca é explorado em
todas as direções
Após algumas gerações 
comunidades envolvidas entre
si  MISTURA  origem a
outras comunidades
SCE-UA
Espaço de busca
Subcomplexo ou subconjunto extraído
aleatoriamente de um complexo 
representa os PAIS
Maior chance
(distribuição triangular
de probabilidades)
aos melhores pais
(menor FO)  geração
de descendentes 
Eles geram outros
complexos
Por vezes, descendentes
são substituídos 
MUTAÇÃO
SCE-UA
Geram-se novos pontos pra
recomeçar o processo até chegar a
um nível satisfatório
SCE-UA
Geram-se novos pontos pra
recomeçar o processo até chegar a
um nível satisfatório
SCE-UA
Geram-se novos pontos pra
recomeçar o processo até chegar a
um nível satisfatório
5
4.5
4
3.5
3
Passo 1
2.5
2
1.5
1
0.5
0
0
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
5
4.5
4
3.5
3
Passo 2
2.5
2
1.5
1
0.5
0
0
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
5
4.5
4
3.5
3
Passo 3
2.5
2
1.5
1
0.5
0
0
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
5
4.5
4
3.5
3
Passo 4
2.5
2
1.5
1
0.5
0
0
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
5
4.5
4
3.5
3
Passo 5
2.5
2
1.5
1
0.5
0
0
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
5
4.5
4
3.5
3
Passo 6
2.5
2
1.5
1
0.5
0
0
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
5
4.5
4
3.5
3
Passo 7
2.5
2
1.5
1
0.5
0
0
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
5
4.5
4
3.5
3
Passo 8
2.5
2
1.5
1
0.5
0
0
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
5
4.5
4
3.5
3
Passo 9
2.5
2
1.5
1
0.5
0
0
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
5
4.5
4
3.5
3
Passo 10
2.5
2
1.5
1
0.5
0
0
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
5
4.5
4
3.5
3
Passo 20
2.5
2
1.5
1
0.5
0
0
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
1 - Geração aleatória de pontos
5
Cada ponto  conj.
de parâmetros
4.5
4
3.5
No caso ao lado 
cada ponto é um par
de parâmetros x1,
x2
3
2.5
2
1.5
Geração aleatória
com uma
distribuição uniforme
1
0.5
0
3:500
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
2 – Hierarquização e formação complexos
5
Hierarquização
segundo as FOs de
cada ponto
4.5
4
Divisão em
complexos 
“comunidades”
3.5
3
2.5
cada um tem
pontos bons e
pontos ruins
2
1.5
Exemplo: complexos
de 4 pontos
1
0.5
0
3:500
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
3 – Evolução dos complexos  sub-complexos
Cada complexo tem
a oportunidade de
evoluir de forma
Independente  as
“comunidades”
evoluem de forma
independente
5
4.5
4
3.5
3
2.5
Nas
“comunidades”
surgem os
“casais”  subcomplexos
2
1.5
1
0.5
0
3:500
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
3 – Evolução dos complexos  sub-complexos
Seleção de “pais”
ou “reprodutores”
 um subcomplexo é
retirado
aleatoriamente de
cada complexo.
5
4.5
4
3.5
3
A probabilidade de
um ponto do
complexo participar
do sub-complexo
é proporcional à
FO (aptidão para
ser “pai”)
2.5
2
1.5
1
0.5
0
3:500
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
3 – Evolução dos complexos  sub-complexos
5
Exemplo:
sub-complexo de 3
pontos extraído de
um complexo de 4
pontos
4.5
4
3.5
3
A probabilidade de
um ponto do
complexo participar
do sub-complexo
é proporcional à
FO (aptidão para
ser “pai”)
2.5
2
1.5
1
0.5
0
3:500
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
3 – Evolução dos complexos  sub-complexos
5
Procedimento de
reprodução 
semelhante ao
método Simplex de
Nelder e Mead
4.5
4
3.5
3
2.5
Define pior
ponto do
sub-complexo
2
1.5
1
0.5
0
3:500
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
3 – Evolução dos complexos  sub-complexos
5
Procedimento de
reprodução 
semelhante ao
método Simplex de
Nelder e Mead
4.5
4
3.5
3
2.5
Define
centróide dos
melhores
pontos
2
1.5
1
0.5
0
3:500
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
3 – Evolução dos complexos  sub-complexos
5
4.5
Passo de reflexão:
distância a = distância b
4
3.5
a
3
Verifica valor da FO
no novo ponto,
se é melhor do que
pior ponto, novo
ponto é aceito,
se não, vai para o
passo de contração.
2.5
b
2
1.5
1
0.5
0
3:500
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
3 – Evolução dos complexos  sub-complexos
5
4.5
Passo de contração:
distância a = distância b
4
a
3.5
b
3
Verifica valor da FO
no novo ponto,
se é melhor do que
pior ponto, novo
ponto é aceito,
se não, cria ponto
aleatório
2.5
2
1.5
1
0.5
0
3:500
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
3 – Evolução dos complexos  sub-complexos
5
4.5
Um novo ponto é gerado
no espaço definido pelos
limites mínimo e máximo
de cada um dos
parâmetros no complexo
4
3.5
3
2.5
2
1.5
1
0.5
0
3:500
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
3 – Evolução dos complexos  sub-complexos
5
Um novo ponto é gerado
no espaço definido pelos
limites mínimo e máximo
de cada um dos
parâmetros no complexo
4.5
4
3.5
Após um n° de
gerações prédeterminado
(β)  comunidades
se encontram
(mistura de
complexos)  nova
divisão é feita
3
2.5
2
1.5
1
0.5
0
3:500
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
4 – Mistura de complexos
5
Mistura de
“comunidades” 
permite o
compartilhamento de
informações que cada
comunidade adquiriu
de forma independente
4.5
4
3.5
3
Complexos
novamente
agrupados 
novo conjunto de
pontos  nova
reorganização...
2.5
2
1.5
1
0.5
0
0
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
4 – Teste de convergência
Pontos da
amostra testados
 caso
satisfaçam um
critério
de convergência
previamente
definido 
algoritmo se
encerra. Caso
contrário 
recomeça no
passo 2.
5
4.5
4
3.5
3
2.5
2
1.5
1
0.5
0
3:500
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
SCE-UA
Pais mais adaptados têm maior probabilidade de gerar filhos
(sobrevivência do mais apto = seleção natural)
Probabilidade de escolha
1
0
Valor da FO
Posição no ranking
1) Classificar os pontos
do complexo em ordem
de FO (ranking)
2) Atribuir probabilidade
de escolha para
participar do subcomplexo segundo a
função do desenho
Exemplo
Sub-Complexo
Complexo
5
5
4.5
4.5
4
4
3.5
3.5
3
3
2.5
2.5
2
2
1.5
1.5
1
1
0.5
0.5
0
0
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
0
0
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
Dois pontos do complexo ficaram fora do sub-complexo. Não
necessariamente os piores pontos ficam fora
Filhos são semelhantes aos pais
Genética: filhos “recebem” cromossomos
dos pais
5
Algoritmo SCE-UA:
No lugar dos “casais” estão
os “sub-complexos”, que são
“casais” de q pontos
4.5
4
3.5
a
3
2.5
b
2
1.5
1
0.5
0
3:50
0
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
Calibração automática com SCE-UA
Cada ponto  valores dos parâmetros escolhidos. FO =
coeficiente de Nash-Sutcliffe. Para ser avaliada  realizar
uma simulação completa (por exemplo, 10 anos de dados
diários).
5
4.5
4
3.5
3
700
2.5
600
calculada
observada
2
Vazão (m3 /s)
500
400
1.5
300
1
200
0.5
100
0
0
0
01/jun/72
01/jul/72
31/jul/72
30/ago/72
29/set/72
29/out/72
28/nov/72
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
Exemplo de aplicação do SCE-UA
COLLISCHONN, W. ; TUCCI, C. E. M. Calibração automática global do
modelo IPH-2. In: Simpósio Brasileiro de Recursos Hídricos, 2001,
Aracaju SE. Anais do Simpósio Brasileiro de Recursos Hídricos, 2001.
Objetivos: testar a robustez, a eficiência e o uso de 2 FOs
diferentes com série contínua de vazões diárias
Robustez  capacidade de encontrar sempre o mesmo
ponto ótimo em uma série de tentativas independentes 
indicativo forte de que o ótimo global foi encontrado
Eficiência capacidade de encontrar o ponto ótimo em um
mínimo número de avaliações da função (iterações
do modelo hidrológico)
Exemplo: teste 1
Como ter certeza de que um método de calibração
pode atingir o ótimo global?
Testá-lo, usando uma série de vazões geradas pelo
próprio modelo  adotando valores de parâmetros
arbitrariemente
O método de calibração tem que atingir estes
valores
I0
Ib
h
Ks
Kbas
Rmax
Alf
50,0
1,0
0,8
5,0
100,0
4,0
2,0
Resultados teste 1
10 aplicações sucessivas  atingido sempre o ótimo
global (conjunto de parâmetros que gerou a série
sintética), em menos do que 10.000 avaliações da FO
140
Valor do parâmetro ao longo do processo
120
Valor do parâmetro I0
Literatura
mostra testes
semelhantes com
métodos
Rosenbrock e
outros, que não
conseguem
superar este
teste
100
I0 = 50
80
60
40
20
0
0
2000
4000
6000
8000
Número de avaliações da função
10000
12000
Teste 2: calibração dados observados
Calibração do modelo IPH-2 (10 vezes)
Calibração
I0
Ib
h
Ks
Kbas
Rmax
Alf
R2
1
36,04
0,46
0,93
7,52
11,11
2,80
19,99
0,85915
2
36,04
0,46
0,93
7,52
11,16
2,80
19,99
0,85915
3
36.03
0.46
0.93
7,52
11,05
2,80
19,99
0,85915
4
35.91
0.46
0,93
7.56
11.95
2.81
19.69
0,85915
5
36.02
0.46
0,93
7.52
11.09
2.80
19.98
0,85915
6
36.04
0.46
0.93
7.52
11.14
2.80
19.99
0,85915
7
36.04
0.46
0,93
7.52
11.12
2.80
19.99
0,85915
8
36.05
0.46
0,93
7.52
11.13
2.80
19.99
0,85915
9
36.03
0.46
0.93
7.52
11.11
2.79
19.99
0,85915
10
36.04
0.46
0.93
7.52
11.16
2.80
19.99
0,85915
Eficaz  indicativo forte de que foi atingido o ótimo global
(apesar de não se conseguir provar)
Otimização multi-objetivo


Objetivo  otimizar, ao mesmo tempo,
várias FOs  vários aspectos da qualidade
do ajuste entre hidrogramas calculados e
observados
Praticamente impossível encontrar um
conjunto de parâmetros que produza, ao
mesmo tempo, valores ótimos de várias
FOs diferentes (imperfeições dos modelos
e dos dados de entrada)
Otimização multi-objetivo


Quase sempre, a opção pela otimização de
uma FO implicará na não otimização de
outra  situação normalmente enfrentada
na calibração de modelos chuva vazão 
bom ajuste dos picos de vazão implica em
maus resultados durante os períodos de
estiagem
A solução não é um ponto, mas uma região
de Pareto, ou região de soluções não
inferiores ou não dominadas
Otimização multi-objetivo
12
11
10
9
F1
F2
Região de Pareto
8
7
F(x)
6
5
4
3
2
1
0
-1
Ótimo F2
Ótimo F1
-2
0
1
2
3
4
5
6
x
Região de
Pareto
7
8
9
10
Otimização multi-objetivo
12
11
10
9
F1
F2
Região de Pareto
8
7
F(x)
6
5
Piorou em F1
4
3
2
1
Mínimo em F2
0
-1
-2
0
1
2
3
4
5
x
x 6= 6
7
8
9
10
Otimização multi-objetivo

A otimização mono-objetivo pode ser
utilizada de forma iterativa para resolver o
problema multi-objetivo:


uso de coeficientes, ou pesos, atribuídos a cada
uma das FOs  a cada calibração  diversos
pontos que definem, aproximadamente, a região de
Pareto (RP) do problema de calibração multiobjetivo
“hierarquização de Pareto” ou Pareto ranking
 encontrar, em apenas um procedimento de
otimização, vários pontos na RP  MOCOM-UA
MOCOM-UA (Multiple-Objective Complex Evolution Universidade do Arizona)




Variação do algoritmo SCE-UA (Duan et al.,
1992), desenvolvido para problemas monoobjetivo
Inicia com a definição dos valores limite dos n
parâmetros
Adota-se o n° de pontos ns (ns conjuntos de
parâmetros)  tamanho da população
São avaliadas nf FOs  matriz de resultados
[ns, nf]
MOCOM-UA
Exemplo de matriz de resultados: n = 1 (1
parâmetro, ns = 3 (três pontos) e nf = 2 (2 FOs)
12
N° do
ponto
1
2
3
11
10
9
F1
F2
Região de Pareto
8
7
F(x)
6
5
4
3
x F1 F2
3
4
7
2
1
0
-1
-2
0
1
2
3
4
5
x
6
7
8
9
10
1
0
9
9
4
1
MOCOM-UA
Exemplo de 2
parâmetros: n = 2
(dois parâmetros,
ns = 5 (cinco
pontos) e nf = 2
(duas FOs)
Ponto Qualquer
(X,Y)
MOCOM-UA


Com os valores das FOs  conjuntos de
parâmetros avaliados e hierarquizados pelos
critérios de dominância e não dominância
Passos:
1. Toma-se os ns pontos ou conjuntos e identifica-se
aqueles que são dominados e os que são não
dominados
2. Aos pontos que são não dominados atribui-se o índice 1
3. Os pontos com o índice 1 são retirados e os pontos
restantes são novamente analisados
4. Aos pontos que são não dominados nesta segunda análise,
atribui-se o índice 2
MOCOM-UA


Com os valores das FOs  conjuntos de
parâmetros avaliados e hierarquizados pelos
critérios de dominância e não dominância
Passos:
5. Os pontos com o índice 2 também são retirados e os
pontos restantes são analisados
6. Os passos se repetem até que se encontre um grupo de
pontos em que não podem ser definidos dominados e não
dominados. Estes pontos recebem o índice MAXR, onde
MAXR é o número de passos necessários para não existir
mais dominância entre os pontos
Quanto menor o índice  mais perto da região
de Pareto
MOCOM-UA: dominância e não dominância
12
11
x = 3 vs x = 7  quem domina quem?
10
9
F1
F2
Região de Pareto
8
7
F(x)
6
5
4
3
2
1
0
-1
-2
0
1
2
3 3
x=
4
5
x
6
x 7= 7
8
9
10
MOCOM-UA: dominância e não dominância
12
11
x = 6 vs x = 7  quem domina quem?
10
9
F1
F2
Região de Pareto
8
7
F(x)
6
5
4
3
2
1
0
-1
-2
0
1
2
3
4
5
x
x =6 6 x 7= 7
8
9
10
MOCOM-UA: evolução dos pontos

Exemplo de 2 parâmetros:



geram-se NMAXR complexos (grupos de
pontos). NMAXR  n° de pontos que recebem
o pior índice na hierarquização
Cada complexo possui n+1 (3 pontos), e a cada
ponto dele é associada uma probabilidade de
que esse ponto participe do processo de
evolução
Cada complexo evolui uma única vez, gerando
NMAXR novos pontos que substituem os
NMAXR piores pontos anteriores, mantendo a
população estável
MOCOM-UA: geração de pontos
Gerados ns =
5 pontos
(conjuntos de
pares X,Y)
MOCOM-UA: hierarquização
Pontos
hierarquizados
(dominância e
não
dominância)
MAXR = 3 e há
somente 1
ponto que
recebe este
índice 
NMAXR = 1
B e C não dominados
B
C
MOCOM-UA: criação de complexos
Gerado
somente 1
complexo
(NMAXR = 1)
São escolhidos
para participar
dele os pontos
B e C (além do
ponto D)
MOCOM-UA: criação de complexos
Gerado
somente 1
complexo
(NMAXR = 1)
MOCOM-UA: evolução dos complexos
Centróide dos
melhores pontos
Método simplex
de Nelder e
Mead
(semelhante ao
algoritmo SCEUA)
MOCOM-UA
reflexão
Ponto R é
aceito se
estiver na
região válida e
se for não
dominado em
relação aos
outros 2 pontos
MOCOM-UA
contração
Caso contrário,
o ponto de
contração K é
aceito
imediatamente
MOCOM-UA
Cada complexo evolui uma única vez, gerando
NMAXR novos pontos que substituem os NMAX
piores pontos anteriores, mantendo a população
estável.
A seguir, pontos novamente analisados e
hierarquizados  processo se repete, até que, na
etapa de hierarquização, todos os pontos recebam o
mesmo índice  não há pontos melhores ou piores.
Nesta situação, normalmente, deverá ter sido
encontrada uma boa amostra de pontos sobre a
região de Pareto
MOCOM-UA
Região de Pareto
Aproximação da região de Pareto
Exemplo IPH 2
Faixa válida dos parâmetros
Parâmetr Unidade
o
Io
mm.Dt-1
Ib
mm.Dt-1
H
Ks
Dt
Ksub
Dt
Rmáx
mm
Alf
-
Valor
mínimo
10
0,1
0,0
0,01
30,0
0,0
0,01
Valor
máximo
300
10
1,0
10,0
40,0
9,0
20,0
2 FO: Erro volume e RMSE
Geração 1
8
7
Erro no volume
6
5
4
3
2
1
0
15000
17000
19000
21000
Soma desvios quadrados
23000
25000
Geração 10
8
7
Erro no volume
6
5
4
3
2
1
0
15000
17000
19000
21000
Soma desvios quadrados
23000
25000
Geração 20
8
7
Erro no volume
6
5
4
3
2
1
0
15000
17000
19000
21000
Soma desvios quadrados
23000
25000
Geração 50
8
7
Erro no volume
6
5
4
3
2
1
0
15000
17000
19000
21000
Soma desvios quadrados
23000
25000
Geração 138
8
7
Erro no volume
6
5
4
3
2
1
0
15000
17000
19000
21000
Soma desvios quadrados
23000
25000
Avaliação da incerteza: usar todos os conjuntos e gerar vários hidrogramas
Propagação da incerteza: Q90 calculada, por exemplo, vai de 8,9 a 10,5
m3.s-1, sendo que a Q90 observada é de 9,1 m3.s-1
Quando se otimiza, aceita-se que
todos os dados são bons, ou seja, a
otimização não faz distinção entre
dados bons e dados ruins
Sugestões de leitura




TUCCI, C. E. M.; COLLISCHONN, W. Ajuste Multiobjetivo dos
Parâmetros de um Modelo Hidrológico. Revista Brasileira de Recursos
Hídrcos, Vol. 8, n. 3, Jul/Set 2003 p. 27-39
BRAVO, J. M.; COLLISCHONN, W; TUCCI, C. E. M. Verificação da
Eficiência e Eficácia de um Algoritmo Evolucionário Multi-objetivo na
Calibração Automática do Modelo Hidrológico IPH II. Revista Brasileira
de Recursos Hídrcos, Vol. 14, n. 3, Jul/Set 2009 p. 37- 50
Yapo, P. O.; Gupta, H. V.; Sorooshian, S. 1998 Multi-objective global
optimization for hydrologic models. Journal of Hydrology, Vol. 204 pp.
83-97.
Sorooshian, S.; Gupta, V. K. 1995 Model calibration In: Singh, V. J.
(editor) Computer models of watershed hydrology. Water Resources
Publications, Highlands Ranch. 1130 p.
Sugestões de leitura


Duan, Q.; Sorooshian, S.; Gupta, V. 1992 Effective and efficient global
optimization for conceptual rainfall-runoff models. Water Resources
Research Vol. 28 No. 4. pp. 1015-1031.
Duan, Q.; Sorooshian, S.; Gupta, V. 1994 Optimal use of the SCE –
UA global optimization method for calibrating watershed models.
Journal of Hydrology, Vol 158 pp. 265-284.
Sugestões de leitura





Bonabeau, E.; Dorigo, M.; Theraulaz, G. 2000
Inspiration for optimization from social insect behaviour.
Nature Vol. 406 July pp.39-42.
Goldberg, D. 1989 Genetic Algorithms in Search,
Optimization and Machine Learning Addison-Wesley, 412
pp.
Klemes, V. 1986 Operational testing of hydrological
simulation models. Hydrological Sciences Journal V. 31 No. 1
pp. 13-24.
Nash e Sutcliffe, 1970 (Journal of Hydrology)
Particle Swarm Optimization
Download

Função objetivo