Vamos montar uma matriz (X) correspondente às leituras dos 10 pluviômetros, mais uma coluna de “uns” referente à constante β₀.

Inicialmente carregar as bibliotecas NumPy, Pandas e Matplotlib.

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt

Abrir o arquivo rainfall.csv, que deve estar no diretório local.

data = pd.read_csv('rainfall.csv', sep=',')

Os dados são armazenados na variável data, do tipo “dataframe”, com a seguinte estrutura:

Para realizar as operações matriciais com os comandos do Numpy, criamos a matriz de dados X extraindo as 10 primeiras colunas da variável data, e a matriz de escoamento Y com a última coluna da variável data.

X = data.values[:,1:11]
Y = data.values[:,11]

A variável X é uma matriz bidimensional de 19 linhas por 10 colunas. Mas a matriz Y é um vetor unidimensional com 19 elementos:

>>> X.shape
(19, 10)
>>> Y.shape
(19,)

Por isso alteramos a dimensionalidade da variável Y com o comando reshape(-1,1):

>>> Y = Y.reshape(-1,1)
>>> Y.shape
(19, 1)

Para o cálculo da estimativa do parâmetro constante (b₀) é necessário inserir uma coluna de “ums” na matriz X e para isso usamos o método ones passando como como argumento o número de linhas da matriz X com o método shape[0]:

>>> ones = np.ones(X.shape[0])
>>> ones
array([1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1.,
       1., 1.])
>>> ones.shape
(19,)
>>> ones = ones.reshape(-1,1)
>>> ones.shape
(19, 1)

A variável ones “aponta” para um vetor unidimensinal de “ums”, por isso usamos o método reshape(-1,1) para transformar em uma matriz com 19 linhas e 1 coluna.

A matriz ones é combinada com a matriz X usando o método hstack([ones, X]):

>>> X = np.hstack([ones, X])

A matriz transposta de X (X_t) é obtida com o método T:

>>> X_t = X.T

O produto matricial de X_t e X é obtido com o método dot(X_t,X) e armazenado na variável X_t_dot_X:

>>> X_t_dot_X = np.dot(X_t,X)

A matriz inversa de X_t_dot_X é calculada com o método linalg.inv(X_t_dot_X) e o resultado armazenado na variável inv_X_t_dot_X:

>>> inv_X_t_dot_X = np.linalg.inv(X_t_dot_X)

O segundo termo da equação é obtido pelo produto matricial das matrizes X_t e Y:

>>> X_t_dot_Y = np.dot(X_t,Y)
>>> X_t_dot_Y = X_t_Y.reshape(-1,1)

E a matriz B com as estimativas dos parâmetros β₀, β₁ ... β₁₀ da equação de regressão pode então ser determinada com o produto matricial:

>>> B = np.dot(inv_X_t_dot_X, X_t_Y)
>>> B
array([[ 7.82347374e+02],
       [ 1.86128582e-01],
       [ 4.83564765e-02],
       [-1.98488741e-02],
       [ 1.91236201e-03],
       [ 1.19622640e-01],
       [ 1.55539685e-01],
       [ 2.31577538e-02],
       [ 1.94764762e-01],
       [ 7.99100377e-02],
       [-4.08511690e-03]])

Cada elemento da matriz B é um dos coeficientes do modelo e portanto o modelo empírico “multivariado” que permite calcular o escoamento anual total da bacia (em mm) a partir dos registros de chuva anual dos 10 pluviômetros pode ser agora expresso pela equação T.2.

Os valores de escoamento da bacia “estimados” pelo modelo são obtidos pelo produto:

>>> Y_est = np.dot(X, B)
>>> Y_est
array([[3264.38629593],
       [3725.6400493 ],
       [3956.259989  ],
       [3889.84426206],
       [3759.56562963],
       [3106.41116674],
       [3068.68275557],
       [2878.20102981],
       [3051.38134815],
       [3884.91244724],
       [3182.54558156],
       [3830.90297952],
       [3324.63158081],
       [3584.66737044],
       [3901.60967752],
       [4625.11081317],
       [3475.56739368],
       [3367.48158331],
       [3941.39804656]])

Os gráficos da figura T.5 mostram os valores de escoamento da bacia estimados pelo modelo em relação aos valores de escoamento observados experimentalmente, os resíduos e a correlação entre valores observados e estimados.

Figura T.5. À esquerda o gráfico do escoamento anual da bacia (azul) e os valores simulados pelo modelo (verde). No centro o gráfico dos resíduos (diferença entre os valores observados e os valores simulados). À direita o gráfico a correlação entre os valores observados experimentalmente e os valores simulados.

Os gráficos da figura T.5 foram gerados com os comandos:

#Exibindo os gráficos

#Ajuste do tamanho
plt.figure(figsize=(15,5))

#Cria subplot com separação horizontal
plt.subplot(1,3,1)

#Extraindo os anos T

T = data.values[:,0]

print("T :", T)

plt.plot(T, Y, 's', c='blue', label='Observado')

plt.plot(T, Y_est, '-', c='green', label='Simulado')

plt.xlabel('Ano')

plt.ylabel('Escoamento da bacia (mm)')

plt.title('Escoamento anual de uma bacia hidrográfica')

plt.legend()

plt.subplot(1,3,2)

#Cálculo da matriz dos resíduos E (erros de previsão)

E = Y - Y_est

plt.plot(T, E, 's', c='red')

plt.hlines(y=0, xmin=T[0], xmax=T[-1], color='blue')

plt.xlabel('Ano')

plt.ylabel('Resíduos (mm)')

plt.title('Resíduos')

plt.subplot(1,3,3)

plt.plot(Y, Y_est, 's', c='blue')

x_lim = plt.xlim()
y_lim = plt.ylim()

print(x_lim, y_lim)

#plt.plot(x_lim, y_lim, c='black')

#Gráfico Identidade
plt.plot( [2600, 4700], [2600, 4700], c='black')

plt.xlabel('Observado')

plt.ylabel('Simulado')

#Aumentar espaço entre os gráficos
plt.subplots_adjust(wspace=0.4)

plt.show()

Water Quality Evaluation of the Yangtze River in China Using Machine Learning Techniques

Visual artificial tongue for identification of various metal ions in mixtures and real water samples: a colorimetric sensor array using off-the-shelf dyes

Medicina e Probabilidade

Vídeos da Acadgild

Principal Component Analysis Tutorial Part 1

Principal Component Analysis Tutorial Part 2

Mod-07 Lec-31 Principal Component Analysis

Mod-07 Lec-32 Regression on Principal Components