lista2-adar-Jestogon
Jéssica Stobienia Gonçalves
2024-10-24
Lista dos exercícios do capítulo 7.14 do ebbok ADAR
Exercício 1
Crie um vetor com os valores de \(e^{x}cosx\) para os valores de \(x=(3,3.1,3.2,...,6)\).
x_1 <- seq(from = 3, to = 6, by = 0.1)
x_1## [1] 3.0 3.1 3.2 3.3 3.4 3.5 3.6 3.7 3.8 3.9 4.0 4.1 4.2 4.3 4.4 4.5 4.6 4.7 4.8
## [20] 4.9 5.0 5.1 5.2 5.3 5.4 5.5 5.6 5.7 5.8 5.9 6.0vetor_1 <- exp(x_1) * cos(x_1)
vetor_1## [1] -19.884531 -22.178753 -24.490697 -26.773182 -28.969238 -31.011186
## [7] -32.819775 -34.303360 -35.357194 -35.862834 -35.687732 -34.685042
## [13] -32.693695 -29.538816 -25.032529 -18.975233 -11.157417 -1.362099
## [19] 10.632038 25.046705 42.099201 61.996630 84.929067 111.061586
## [25] 140.525075 173.405776 209.733494 249.468441 292.486707 338.564378
## [31] 387.360340Exercício 2
Crie os seguintes vetores:
- \((0.1³×0.2¹,0.1⁶×0.2⁴,...,0.1^{36}×0.2^{34})\)
(x_2 <- seq(from = 3, to = 36, by = 3))## [1] 3 6 9 12 15 18 21 24 27 30 33 36(y_2 <- x_2 - 2)## [1] 1 4 7 10 13 16 19 22 25 28 31 34vetor_2 <- 0.1^(x_2) * 0.2^(y_2)
vetor_2## [1] 2.000000e-04 1.600000e-09 1.280000e-14 1.024000e-19 8.192000e-25
## [6] 6.553600e-30 5.242880e-35 4.194304e-40 3.355443e-45 2.684355e-50
## [11] 2.147484e-55 1.717987e-60- \((2,{2^{2}}/2,{2^{3}}/3,...,{2^{25}}/25)\)
(x_3 <- seq(from = 2, to = 25, by = 1))## [1] 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25vetor_3 <- 2^(x_3) / x_3
vetor_3## [1] 2.000000e+00 2.666667e+00 4.000000e+00 6.400000e+00 1.066667e+01
## [6] 1.828571e+01 3.200000e+01 5.688889e+01 1.024000e+02 1.861818e+02
## [11] 3.413333e+02 6.301538e+02 1.170286e+03 2.184533e+03 4.096000e+03
## [16] 7.710118e+03 1.456356e+04 2.759411e+04 5.242880e+04 9.986438e+04
## [21] 1.906502e+05 3.647221e+05 6.990507e+05 1.342177e+06Exercício 3
Reproduza a criação do vetor dias da semana \((dds)\) mostrado abaixo:
dias <- c("domingo", "segunda", "terca", "quarta", "quinta", "sexta", "sabado")
dds <- c(1:7)
names(dds) <- dias
dds## domingo segunda terca quarta quinta sexta sabado
## 1 2 3 4 5 6 7Exercício 4
Interprete o resultado da seguinte operação:
dds_nums <- c(
5L, 2L, 5L, 2L, 7L, 7L, 2L,
6L, 6L, 3L, 7L, 1L, 2L, 2L,
5L, 7L, 3L, 2L, 6L, 4L, 1L
)
names(dds)[dds_nums]## [1] "quinta" "segunda" "quinta" "segunda" "sabado" "sabado" "segunda"
## [8] "sexta" "sexta" "terca" "sabado" "domingo" "segunda" "segunda"
## [15] "quinta" "sabado" "terca" "segunda" "sexta" "quarta" "domingo"Interpretação: através da indexação os números do dds_nums representam os dias da semana. _______________________________________________________________________________ ## Exercício 5
- Escreva o código necessário para determinar o vetor lógico indicando quais números são pares na sequência de valores 85, 79, 70, 6, 32, 8, 17, 93, 81, 76.
x_4 <- c(85, 79, 70, 6, 32, 8, 17, 93, 81, 76)
pares <- x_4 %% 2 == 0
pares## [1] FALSE FALSE TRUE TRUE TRUE TRUE FALSE FALSE FALSE TRUE- Calcule o total de números ímpares.
x_4 <- c(85, 79, 70, 6, 32, 8, 17, 93, 81, 76)
impares <- x_4 %% 2 != 0
total_impares <- sum(impares)
total_impares## [1] 5Exercício 6
Para um ano ser bissexto, ele deve ser:
divísivel por 4 (a divisão é exata com resto igual a zero);
não pode ser divisível por 100 (a divisão não é exata, ou seja, o resto é diferente de zero);
pode ser que seja divisível por 400: caso seja divisível por 400, a divisão deve ser exata, deixando o resto igual a zero.
Com os critérios definidos acima, construa o código para:
- Verificar se os anos 1894, 1947, 1901, 1992, 1925, 2014, 1993, 1996, 1984, 1897, 2100, 2300 são bissextos.
anos <- c(1894, 1947, 1901, 1992, 1925, 2014, 1993, 1996, 1984, 1897, 2100, 2300)
bissextos <- (anos %% 4 == 0 & anos %% 100 != 0) | (anos %% 400 == 0)
bissextos## [1] FALSE FALSE FALSE TRUE FALSE FALSE FALSE TRUE TRUE FALSE FALSE FALSE- Mostre quais anos são bissextos.
anos[bissextos]## [1] 1992 1996 1984- Usando o código para verificar se o ano é bissexto, gere um vetor nomeado ndias com o número de dias do ano para os anos do item (a).
ndias <- ifelse(test = bissextos, yes = 366, no = 365)
names(ndias) <- anos
ndias## 1894 1947 1901 1992 1925 2014 1993 1996 1984 1897 2100 2300
## 365 365 365 366 365 365 365 366 366 365 365 365- Programe como obter o total de anos com 366 dias?
total_bissextos <- sum(bissextos)
total_bissextos## [1] 3Exercício 7
Quais códigos para gerar os seguintes dados:
## [1] -30 0 5 10 15 20 30x_5 <- c(-30, seq(0, 20, by = 5), 30)
x_5## [1] -30 0 5 10 15 20 30## [1] 1.0 0.8 0.6 0.4 0.2 0.0x_6 <- seq(1, 0, by = -0.2)
x_6## [1] 1.0 0.8 0.6 0.4 0.2 0.0## [1] -3.1415927 -2.4434610 -1.7453293 -1.0471976 -0.3490659 0.3490659
## [7] 1.0471976 1.7453293 2.4434610 3.1415927x_7 <- seq(-pi, pi, by = 0.6981317)
x_7## [1] -3.1415927 -2.4434610 -1.7453293 -1.0471976 -0.3490659 0.3490658
## [7] 1.0471975 1.7453292 2.4434609 3.1415926## [1] -1 -1 0 0 0 1 1 1 1 2 2 2 2 2 3 3 3 3 3 3 4 4 4 4 4
## [26] 4 4 5 5 5 5 5 5 5 5x_8 <- rep(-1:5, times = 2:8)
x_8## [1] -1 -1 0 0 0 1 1 1 1 2 2 2 2 2 3 3 3 3 3 3 4 4 4 4 4
## [26] 4 4 5 5 5 5 5 5 5 5## [1] 5 5 5 5 5 4 4 4 4 3 3 3 2 2 1 2 2 3 3 3 4 4 4 4 5 5 5 5 5x_9 <- c(rep(5:2, times = 5:2), 1, rep(2:5, times = 2:5))
x_9## [1] 5 5 5 5 5 4 4 4 4 3 3 3 2 2 1 2 2 3 3 3 4 4 4 4 5 5 5 5 5Exercício 8
Usando o mesmo código para solução em todos os itens abaixo, obtenha as seguintes sequências usando os vetores fornecidos.
- \(v3 = (11, 0.25, 7, 2)\)
## [1] 1 2 3 4v3 <- c(11, 0.25, 7, 2)
seq_v3 <- seq_along(v3)
seq_v3## [1] 1 2 3 4- \(v2 = (11, 0.25)\)
## [1] 1 2v2 <- c(11, 0.25)
seq_v2 <- seq_along(v2)
seq_v2## [1] 1 2- \(v1 = (11)\)
## [1] 1v1 <- c(11)
seq_v1 <- seq_along(v1)
seq_v1## [1] 1- \(v0 = ()\)
## integer(0)v0 <- c()
seq_v0 <- seq_along(v0)
seq_v0## integer(0)Exercício 9
Considere os seguintes dados horários de temperatura do ar \((T_{ar})\) registrados em duas estações meteorológicas, entre as 0 e 23 horas de um dado dia.
tar_est1 <- c(
14.92, 14.61, 14.32, 14.07, 13.84, 13.65, 13.56, 13.97, 15.08,
16.5, 17.88, 19.08, 20.02, 20.66, 21.01, 21.05, 20.76, 20.05,
18.77, 17.51, 16.67, 16.11, 15.66, 15.27
)
tar_est2 <- c(
13.13, 13.01, 12.93, 12.87, 12.82, 12.81, 13.2, 14.22, 15.77,
17.49, 19.2, 20.57, 21.49, 22.01, 22.03, 21.71, 20.84, 18.77,
16.54, 15.13, 14.34, 13.81, 13.49, 13.28
)Dica: Faça um gráfico para visualizar as temperaturas das duas estações. Isso facilitará a solução.- Determine a média diária da \((T_{ar})\) das duas estações arrendondando para uma casa decimal. Salve a média de cada estação nas variáveis tmed_est1 e tmed_est2.
(tmed_est1 <- round(mean(tar_est1), 1))## [1] 16.9(tmed_est2 <- round(mean(tar_est2), 1))## [1] 16.3- Utilizando as variáveis do item anterior, verifique usando comparação lógica, em qual estação o ambiente é mais quente?
medias <- c(tmed_est1, tmed_est2)
names(medias) <- c("tmed_est1", "tmed_est2")
medias[which.max(medias)]## tmed_est1
## 16.9- Obtenha a 3ª maior temperatura do dia em cada estação.
(terceira_maiorT_est1 <- sort(tar_est1, decreasing = TRUE)[3])## [1] 20.76(terceira_maiorT_est2 <- sort(tar_est2, decreasing = TRUE)[3])## [1] 21.71- Calcule a amplitude térmica diária (\(ATD=Tmax−Tmin\), onde \(Tmax\): temperatura máxima dária e \(Tmin\): temperatura mínima dária) das estações 1 e 2, arrendondando para uma casa decimal. Salve os resultados nas variaveis atd_est1 e atd_est2.
Tmax2 <- max(tar_est2)
Tmax1 <- max(tar_est1)
Tmin2 <- min(tar_est2)
Tmin1 <- min(tar_est1)
(atd_est1 <- round(Tmax1 - Tmin1, 1))## [1] 7.5(atd_est2 <- round(Tmax2 - Tmin2, 1))## [1] 9.2- Qual o horário de ocorrência das temperaturas máximas e mínimas em cada estação? Salve os resultados nas variáveis hmax_est{i} e hmin_est{i}, com \(i=1,2\).
temps <- list(tar_est1, tar_est2)
hmax_est <- list()
hmin_est <- list()
for (i in 1:2) {
hmax_est[[i]] <- which.max(temps[[i]])
hmin_est[[i]] <- which.min(temps[[i]])
}
names(hmax_est) <- paste0("hmax_ext", 1:2)
names(hmin_est) <- paste0("hmim_ext", 1:2)
hmax_est## $hmax_ext1
## [1] 16
##
## $hmax_ext2
## [1] 15hmin_est## $hmim_ext1
## [1] 7
##
## $hmim_ext2
## [1] 6#salvar as variáveis separadamente:
for (i in 1:2) {
assign(paste0("hmax_est", i), which.max(temps[[i]]))
assign(paste0("hmin_est", i), which.min(temps[[i]]))
}- Quando tar_est2 é maior que tar_est1 qual a maior diferença absoluta de temperatura entre as duas estações?
diferencas <- abs(tar_est2 - tar_est1)
dif_valida <- diferencas[tar_est2 > tar_est1]
maior_dif <- max(dif_valida)
maior_dif## [1] 1.49- Qual a hora correspondende a ocorrência daquela maior diferença absoluta de temperatura obtida no item anterior?
horas <- 0:23
indice_maior_dif <- which(diferencas == maior_dif & tar_est2 > tar_est1)
hora_maior_dif <- horas[indice_maior_dif]
hora_maior_dif## [1] 11maior_dif## [1] 1.49- O horário do pôr do sol pode ser estimado a partir da \(T_{ar}\). Para o período após o horário de ocorrência da \(T_{max}\) determina-se em qual hora ocorre a maior queda de \(T_{ar}\) em relação a hora anterior. Estime o horário do pôr do sol para as duas estações (hps_est{i}).
horario_ps <- list()
for (i in 1:2){
queda_temp <- diff(temps[[i]][(hmax_est[[i]] + 1):length(temps[[i]])])
maior_queda <- which.min(queda_temp)
horario_ps[[i]] <- (hmax_est[[i]] + maior_queda)
}
names(horario_ps) <- paste0("horario_ps", 1:2)
horario_ps## $horario_ps1
## [1] 18
##
## $horario_ps2
## [1] 18- Em quais horas as temperaturas das duas estações estão mais próximas uma da outra, ou seja com menos de 0,5°C de diferença ?
(horas_proximas <- horas[diferencas < 0.5])## [1] 6 7 16- Calcule a temperatura média diária usando os seguintes métodos para estação 2. \(T_{med1}=(T_{max}+T_{min})/2\), salvando em uma variável tar_met1. \(T_{med2}=(T_{max}+T_{min}+T_{9}+2T_{21})/5\), salvando em uma variável tar_met2. \(T_{med3}=(T_{7}+T_{14}+2T_{21})/4\), salvando em uma variável tar_met2.
T9 <- tar_est2[9]
T21 <- tar_est2[21]
T7 <- tar_est2[7]
T14 <- tar_est2[14]
tar_met1 <- (Tmax2 + Tmin2) / 2
tar_met2 <- (Tmax2 + Tmin2 + T9 + 2 * T21) / 5
tar_met3 <- (T7 + T14 + 2 * T21) / 4
tar_met1## [1] 17.42tar_met2## [1] 15.858tar_met3## [1] 15.9725- Compare este resultados com aqueles obtidos no item a. Qual é melhor?
valor_dif_medias2 <- abs(c(tar_met1, tar_met2, tar_met3) - tmed_est2)
names(valor_dif_medias2) <- c("tar_met1", "tar_met2", "tar_met3")
barplot(valor_dif_medias2,
xlab = "Médias alternativas",
ylab = "Diferença de Temperatura com a média padrão ºC",
main = "Comparação entre as médias calculadas para a Estação 2",
ylim = c(0, max(valor_dif_medias2) + 1))
valor_dif_medias2[which.min(valor_dif_medias2)]## tar_met3
## 0.3275Exercício 10
Calcule a temperatura do ar horária usando o modelo de onda para representação do ciclo diário da temperatura do ar, descrito pelas seguintes equações:
\[\left\{\begin{matrix} h_{Tmin}\leqslant h < h_{Tmax}, \ T_{calc}= \bar{T} - Acos(arg) \\ demais \ horas,\ T_{calc}= \bar{T} + Acos(arg) \end{matrix}\right.\]
onde:
\[\bar{T}=\left ( \frac{T_{max}+T_{min}}{2} \right )\]
e
\[A=\left ( \frac{T_{max}-T_{min}}{2} \right )\]
O argumento do cosseno \((arg)\) é definido por:
\[\left\{\begin{matrix} h < h_{Tmin}, \ arg=\left ( \frac{h+10}{10+h_{Tmin}} \right ) \\ h_{Tmin}\leqslant h < h_{Tmax},\ arg=\left ( \frac{h-h_{Tmin}}{14-h_{Tmin}} \right ) \\ h > h_{Tmax}, \ arg=\left ( \frac{14-h}{10+h_{Tmin}} \right ) \end{matrix}\right.\]
- Aplique o método acima para estação 1 e 2 substituindo os valores de \(Tmax, Tmin, h_{Tmax}, h_{Tmin}\),
calcula_Tcalc <- function(Tmax, Tmin, hmax_est, hmin_est) {
T_media <- (Tmax + Tmin) / 2
A <- (Tmax - Tmin) / 2
T_calc <- numeric(24)
for (h in 0:23) {
if (h < hmin_est) {
arg <- (h + 10) / (10 + hmin_est)
} else if (h >= hmin_est & h < hmax_est) {
arg <- (h - hmin_est) / (14 - hmin_est)
} else if (h > hmax_est){
arg <- (14 - h) / (10 + hmin_est)
}
if (h >= hmin_est & h < hmax_est) {
T_calc[h + 1] <- T_media - A * cos(pi * arg)
} else {
T_calc[h + 1] <- T_media + A * cos(pi * arg)
}
}
return(T_calc)
}
# Aplicando o método para as duas estações
T_calc1 <- calcula_Tcalc(Tmax1, Tmin1, hmax_est1, hmin_est1)
T_calc1## [1] 16.28013 15.63571 15.04813 14.53741 14.12094 13.81289 13.62377 13.56000
## [9] 13.93087 14.97003 16.47166 18.13834 19.63997 20.67913 21.05000 20.67913
## [17] 13.93087 20.48906 20.07259 19.56187 18.97429 18.32987 17.65055 16.95945T_calc2 <- calcula_Tcalc(Tmax2, Tmin2, hmax_est2, hmin_est2)
T_calc2## [1] 15.65583 14.85882 14.16024 13.58693 13.16092 12.89858 12.81000 13.16092
## [9] 14.16024 15.65583 17.42000 19.18417 20.67976 21.67908 22.03000 12.81000
## [17] 21.67908 21.25307 20.67976 19.98118 19.18417 18.31937 17.42000 16.52063# Plotando T_calc1 e T_calc2
plot(horas, T_calc1,
type = "o",
col = "blue",
pch = 16,
ylim = range(c(T_calc1, T_calc2)),
xlab = "Hora",
ylab = "Temperatura Calculada (°C)",
main = "Temperatura Horária Calculada para Estações 1 e 2",
xaxt = "n")
axis(1,
at = horas,
labels = horas)
lines(horas, T_calc2,
type = "o",
col = "red",
pch = 16)
legend("topright",
legend = c("T_calc1 (Estação 1)", "T_calc2 (Estação 2)"),
col = c("blue", "red"),
pch = 16,
lty = 1)
- Calcule o RMSE nos dois casos. \(RMSE=\sqrt{\frac{1}{n}\sum_{i=1}^{n}\left ( T_{calc}-T_{obs} \right )^{2}}\)
rmse <- function(T_obs, T_calc) {
sqrt(mean((T_calc - T_obs) ^ 2))
}
RMSE1 <- rmse(tar_est1, T_calc1)
RMSE1## [1] 1.83183RMSE2 <- rmse(tar_est2, T_calc2)
RMSE2## [1] 3.038173- Calcule a correlação \((r)\) nos dois casos. A barra representa a média aritmética. Confira seu resultado com a saída da função cor(tar_obs, tar_calc).
\[r=\frac{\sum_{i=1}^{n}\left ( T_{obs}-\bar{T}_{obs} \right )\left ( T_{calc}-\bar{T}_{calc} \right )}{\sqrt{\sum_{i=1}^{n}\left ( T_{obs}-\bar{T}_{obs} \right )^{2}\left ( T_{calc}-\bar{T}_{calc} \right )^{2}}}\]
calcula_correlacao <- function(tar_est, tmed_est, T_calc) {
mean_calc <- mean(T_calc)
numerador <- sum((tar_est - tmed_est) * (T_calc - mean_calc))
denominador <- sqrt(sum((tar_est - tmed_est)^2) * sum((T_calc - mean_calc)^2))
# Calculando a correlação
r_manual <- numerador / denominador
# Retornando o valor da correlação
return(r_manual)
}
(corr1_R <- cor(tar_est1, T_calc1))## [1] 0.7552633(correlacao1 <- calcula_correlacao(tar_est1, tmed_est1, T_calc1))## [1] 0.755231(corr2_R <- cor(tar_est2, T_calc2))## [1] 0.6144771(correlacao2 <- calcula_correlacao(tar_est2, tmed_est2, T_calc2))## [1] 0.6144741dif_correlacoes <- abs(c(correlacao1, corr1_R, correlacao2, corr2_R))
names(dif_correlacoes) <- c("correlacao1", "corr1_R", "correlacao2", "corr2_R")
barplot(dif_correlacoes,
xlab = "Correlação calculada e função correlação do R",
main = "Comparação entre correlações",
ylim = c(0, max(dif_correlacoes) + 0.2))
_______________________________________________________________________________
## Exercício 11
Os dados abaixo são de precipitação horária de um evento severo ocorrido em 03/12/2012 em Santa Maria-RS.
| hora | prec |
|---|---|
| 9 | 0.0 |
| 10 | 0.0 |
| 11 | 0.0 |
| 12 | 0.0 |
| 13 | 0.0 |
| 14 | 0.0 |
| 15 | 0.0 |
| 16 | 21.4 |
| 17 | 41.2 |
| 18 | 2.6 |
| 19 | 1.0 |
| 20 | 0.4 |
| 21 | 0.0 |
- Como seria o código para determinar a soma cumulativa da precipitação horária? Salve o resultado em um vetor chamado prec_acum. Interprete o resultado de c(NA, diff(prec_acum)).
prec_acum <- cumsum(prec)
prec_acum## [1] 0.0 0.0 0.0 0.0 0.0 0.0 0.0 21.4 62.6 65.2 66.2 66.6 66.6Interpretação de c(NA, diff(prec_acum)): a função diff calcula a diferença entre elementos consecutivos no vetor prec_acum a expressão acima está criando um vetor: adicionando NA no primeiro elemento e mostrando a diferença entre os elementos do vetor prec_acum, ou seja, desfazendo a função cumsum.
(c(NA, diff(prec_acum)))## [1] NA 0.0 0.0 0.0 0.0 0.0 0.0 21.4 41.2 2.6 1.0 0.4 0.0- Mostre o código para encontrar o horário de ocorrência da precipitação máxima?
indice_max_prec <- which.max(prec)
hora_max_prec <- hora[indice_max_prec]
hora_max_prec## [1] 17- Mostre o código para obter a hora de início e fim do evento de precipitação severa. Qual foi a duração do evento? OBS: De acordo com o manual do observador meteorológico (INMET) Chuva de 1.1 a 5mm/h é classificada como FRACA, de 5.1 à 60mm/h MODERADA e acima de 60mm/h FORTE. Utilizando a classificação de intensidade de precipitação MODERADA segue a resolução do exercío…
indices_moderada <- which(prec >= 5.1 & prec <= 60)
inicio_evento <- hora[min(indices_moderada)]
fim_evento <- hora[max(indices_moderada)]
duracao_evento <- fim_evento - inicio_evento
inicio_evento## [1] 16fim_evento## [1] 17duracao_evento## [1] 1- Qual foi a precipitação total do evento? Quanto da precipitação total do evento, em %, ocorreu até às 17 h?
indices_evento <- which(prec > 0)
prec_total_evento <- sum(prec[indices_evento])
prec_ate_17h <- sum(prec[hora <= 17])
percentual_ate_17h <- (prec_ate_17h / prec_total_evento) * 100
prec_total_evento## [1] 66.6prec_ate_17h## [1] 62.6percentual_ate_17h## [1] 93.99399Exercício 12
Considere o vetor x definido pelos números descritos abaixo. Mostre como encontrar o primeiro número positivo localizado após o último número negativo. Por exemplo, seja o vetor z definido pelos valores (11, 10, 15, 2, 6, -15, -10, -22, -8, 5, 7, 2, 12, 8, 4, 1, 3, -3, -1, 30, 14). Os valores selecionados seriam 5 e 30.
z <- c(11, 10, 15, 2, 6, -15, -10, -22, -8, 5, 7, 2, 12, 8, 4, 1, 3, -3, -1, 30, 14)
z_ordenado <- sort(z)
indice_ultimo_negativo <- max(which(z_ordenado < 0))
ultimo_negativo <- z_ordenado[indice_ultimo_negativo]
indice_primeiro_positivo <- min(which(z_ordenado > 0))
primeiro_positivo <- z_ordenado[indice_primeiro_positivo]
z_ordenado## [1] -22 -15 -10 -8 -3 -1 1 2 2 3 4 5 6 7 8 10 11 12 14
## [20] 15 30ultimo_negativo## [1] -1primeiro_positivo## [1] 1Exercício 13
Considerando o vetor prec com valores de precipitação diária indicado abaixo. Escreva o código para resolver as seguintes tarefas.
prec_13 <- c(
0, 0, 0, 0.8, 0, 0.01, 0.75, 0,
0, 0, 0, 0.35, 0.08, 0, 0, 0, 0, 0.31, 0, 3.57, 12.17, 0, 0,
0, 0.04, 3.16, 0, 0.95, 0.79, 0, 0, 0, 0, 0, 3.51, 0, 0, 0.16,
0, 0, 8.16, 0.54, 4.39, 1.24, 0, 0, 0, 0, 0, 2.43, 0, 0, 0, 0,
0, 7.18, 0, 0, 0.26, 0, 0, 0.28, 0, 0, 0.09, 0.38, 0, 0, 0, 0,
0, 0, 0.51, 0, 0, 0, 0, 0, 0, 0.67, 0, 0, 0, 0, 0.15, 0, 0.82,
0, 0, 0, 0, 0, 0, 0, 0, 0.37, 0, 0.58, 4.95, 0, 0, 0, 0, 0, 7.68,
0, 0, 0.37, 0, 1.56, 0, 0, 0, 0.34, 0.48, 0, 4.21, 2.28, 4.3,
0, 3.38, 0, 0, 0, 0, 7.28, 0, 4.89, 3.91, 0, 0, 0, 0, 0, 0, 2.93,
0, 2.49, 0.77, 0, 2.9, 3.53, 0.83, 0, 0, 0, 0.94, 0.59, 0, 0,
0, 0, 0.04, 0, 0.65, 0, 0, 0, 6.23, 0.09, 0, 0.66, 0, 0, 0, 4.42,
0, 0, 0, 0.84, 0, 0, 0, 0, 0, 0.09, 0, 0, 0.08, 0, 0.66, 0, 0,
0, 0.06, 0, 0, 0, 3.28, 0, 0.8, 5.69, 0.8, 0
)- Quantos dias ocorreram no intervalo 0 < prec_13 < 0.25?
dias_intervalo <- sum(prec_13 > 0 & prec_13 < 0.25)
dias_intervalo## [1] 11- Substitua os valores de chuva registrados no intervalo 0 < prec_13 < 0.25 por 0.
prec_13[prec_13 > 0 & prec_13 < 0.25] <- 0
prec_13## [1] 0.00 0.00 0.00 0.80 0.00 0.00 0.75 0.00 0.00 0.00 0.00 0.35
## [13] 0.00 0.00 0.00 0.00 0.00 0.31 0.00 3.57 12.17 0.00 0.00 0.00
## [25] 0.00 3.16 0.00 0.95 0.79 0.00 0.00 0.00 0.00 0.00 3.51 0.00
## [37] 0.00 0.00 0.00 0.00 8.16 0.54 4.39 1.24 0.00 0.00 0.00 0.00
## [49] 0.00 2.43 0.00 0.00 0.00 0.00 0.00 7.18 0.00 0.00 0.26 0.00
## [61] 0.00 0.28 0.00 0.00 0.00 0.38 0.00 0.00 0.00 0.00 0.00 0.00
## [73] 0.51 0.00 0.00 0.00 0.00 0.00 0.00 0.67 0.00 0.00 0.00 0.00
## [85] 0.00 0.00 0.82 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.37
## [97] 0.00 0.58 4.95 0.00 0.00 0.00 0.00 0.00 7.68 0.00 0.00 0.37
## [109] 0.00 1.56 0.00 0.00 0.00 0.34 0.48 0.00 4.21 2.28 4.30 0.00
## [121] 3.38 0.00 0.00 0.00 0.00 7.28 0.00 4.89 3.91 0.00 0.00 0.00
## [133] 0.00 0.00 0.00 2.93 0.00 2.49 0.77 0.00 2.90 3.53 0.83 0.00
## [145] 0.00 0.00 0.94 0.59 0.00 0.00 0.00 0.00 0.00 0.00 0.65 0.00
## [157] 0.00 0.00 6.23 0.00 0.00 0.66 0.00 0.00 0.00 4.42 0.00 0.00
## [169] 0.00 0.84 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00
## [181] 0.66 0.00 0.00 0.00 0.00 0.00 0.00 0.00 3.28 0.00 0.80 5.69
## [193] 0.80 0.00- Crie um vetor denominado prec01 indicando o estado da precipitação (chuvoso = 1, seco = 0) baseado no limiar de 0.25 mm para detecção de chuva pelo pluviômetro.
prec01 <- ifelse(prec_13 >= 0.25, 1, 0)
prec01## [1] 0 0 0 1 0 0 1 0 0 0 0 1 0 0 0 0 0 1 0 1 1 0 0 0 0 1 0 1 1 0 0 0 0 0 1 0 0
## [38] 0 0 0 1 1 1 1 0 0 0 0 0 1 0 0 0 0 0 1 0 0 1 0 0 1 0 0 0 1 0 0 0 0 0 0 1 0
## [75] 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 1 0 1 1 0 0 0 0 0 1 0 0 1 0 1 0
## [112] 0 0 1 1 0 1 1 1 0 1 0 0 0 0 1 0 1 1 0 0 0 0 0 0 1 0 1 1 0 1 1 1 0 0 0 1 1
## [149] 0 0 0 0 0 0 1 0 0 0 1 0 0 1 0 0 0 1 0 0 0 1 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0
## [186] 0 0 0 1 0 1 1 1 0- Qual a probabilidade de chuva dessa série de precipitação diária?
dias_chuvosos <- sum(prec_13 >= 0.25)
total_dias <- length(prec_13)
probabilidade_chuva <- dias_chuvosos / total_dias
probabilidade_chuva## [1] 0.2835052- Qual a probabilidade de chover dois dias consecutivos (p11)? Calcule a probabilidade de chover em qualquer um de dois dias consecutivos (p01 + p10)?
p11 <- sum(prec01[-length(prec01)] == 1 & prec01[-1] == 1) / (length(prec01) - 1)
p01 <- sum(prec01[-length(prec01)] == 0 & prec01[-1] == 1) / (length(prec01) - 1)
p10 <- sum(prec01[-length(prec01)] == 1 & prec01[-1] == 0) / (length(prec01) - 1)
p01_p10 <- p01 + p10
p11## [1] 0.08290155p01_p10## [1] 0.4041451- Determine a duração de cada evento chuvoso (número de dias consecutivos).
duracoes <- rle(prec01)
duracao_eventos <- duracoes$lengths[duracoes$values == 1]
names(duracao_eventos) <- paste0("evento", seq_along(duracao_eventos))
duracao_eventos## evento1 evento2 evento3 evento4 evento5 evento6 evento7 evento8
## 1 1 1 1 2 1 2 1
## evento9 evento10 evento11 evento12 evento13 evento14 evento15 evento16
## 4 1 1 1 1 1 1 1
## evento17 evento18 evento19 evento20 evento21 evento22 evento23 evento24
## 1 1 2 1 1 1 2 3
## evento25 evento26 evento27 evento28 evento29 evento30 evento31 evento32
## 1 1 2 1 2 3 2 1
## evento33 evento34 evento35 evento36 evento37 evento38 evento39
## 1 1 1 1 1 1 3Exercício 14
dados <- c(
NA, NA, 27L, 7L, 4L, 0L, 26L, 15L, 25L, NA,
NA, NA, NA, 6L, 29L, 18L, 17L, 23L, 20L, 1L,
30L, 13L, NA, NA, NA, NA, NA, NA, NA, 19L
)- Como você pode codificar a obtenção de um vetor com zeros nos valores válidos e com números sequenciais dentro das falhas?
resultado <- ifelse(is.na(dados), seq_along(dados[dados]), 0)
resultado## [1] 1 2 0 0 0 0 0 0 0 10 11 12 13 0 0 0 0 0 0 0 0 0 23 24 25
## [26] 26 27 28 29 0- Como a partir do vetor resultante do item anterior você pode obter um vetor cujo os valores dentro das falhas indique a ordem de ocorrência da falha.
numero_na <- seq(1:sum(is.na(dados)))
resultado.na <- resultado[resultado != 0]
names(resultado.na) <- numero_na
resultado.na## 1 2 3 4 5 6 7 8 9 10 11 12 13
## 1 2 10 11 12 13 23 24 25 26 27 28 29- Qual o código para determinar o tamanho de cada falha?
tamanho_falhas <- rle(is.na(dados))$lengths[rle(is.na(dados))$values]
tamanho_falhas## [1] 2 4 7- Como determinar o tamanho da maior falha?
maior_tamanho_falhas <- max(tamanho_falhas)
maior_tamanho_falhas## [1] 7Exercício 15
Para os valores de velocidade \((V_h)\) e a direção do vento \((θ)\) (na convenção meteorológica): \[V_h = (10, 10, 10, 10, 14.142, 14.142, 14.142, 14.142, 0)\] \[θ = (270, 180, 360, 90, 225, 315, 135, 45, 0)\] a. Determine as componentes zonal e meridional do vento. \[u=−V_h⋅sin(θ_{rad})\] \[v=−V_h⋅cos(θ_{rad})\]
V_h <- c(10, 10, 10, 10, 14.142, 14.142, 14.142, 14.142, 0)
theta <- c(270, 180, 360, 90, 225, 315, 135, 45, 0)
theta_rad <- theta * (pi / 180)
u <- -V_h * sin(theta_rad)
v <- -V_h * cos(theta_rad)
tabela_u_v <- data.frame(V_h, theta, u, v)
tabela_u_v## V_h theta u v
## 1 10.000 270 1.000000e+01 1.836970e-15
## 2 10.000 180 -1.224647e-15 1.000000e+01
## 3 10.000 360 2.449294e-15 -1.000000e+01
## 4 10.000 90 -1.000000e+01 -6.123234e-16
## 5 14.142 225 9.999904e+00 9.999904e+00
## 6 14.142 315 9.999904e+00 -9.999904e+00
## 7 14.142 135 -9.999904e+00 9.999904e+00
## 8 14.142 45 -9.999904e+00 -9.999904e+00
## 9 0.000 0 0.000000e+00 0.000000e+00- Faça os cálculos necessários para reconstruir \(V_h\) e \(θ\), a partir de \(u\) e \(v\) determinados no item a. Por convenção, a direção do vento \(θ\) em condições calmas \((V_h < 0.5 ms^{-1})\) é assumida como \(0°\). \[V_h=\sqrt{(u^2+v^2)}\] \[θ_{mat}=atan2(−u,−v)⋅\frac{180}{π}\] θ={θmat+360seθmat<00seu=0,v=0ouVh<0.5 \[θ = \left\{\begin{matrix} θ_{mat}+360 \ se \ θ_{mat} < 0 \\ 0 \ se \ u=0, \ v=0 \ ou \ V_h<0.5 \end{matrix}\right.\]
A tabela abaixo apresenta o resultado esperado para as variáveis derivadas.
| u | v | w_s | w_d | wd_uv | dir |
|---|---|---|---|---|---|
| 10 | 0 | 10.000 | 270 | 270 | Oeste |
| 0 | 10 | 10.000 | 180 | 180 | Sul |
| 0 | -10 | 10.000 | 360 | 360 | Norte |
| -10 | 0 | 10.000 | 90 | 90 | Leste |
| 10 | 10 | 14.142 | 225 | 225 | Sudoeste |
| 10 | -10 | 14.142 | 315 | 315 | Noroeste |
| -10 | 10 | 14.142 | 135 | 135 | Sudeste |
| -10 | -10 | 14.142 | 45 | 45 | Nordeste |
| 0 | 0 | 0.000 | 0 | 0 | Calmo |
Dica: ver figura abaixo.
V_h2 <- sqrt(u^2 + v^2)
theta_mat <- atan2(-u, -v) * (180 / pi)
theta2 <- ifelse(theta_mat < 0, theta_mat + 360, theta_mat)
theta2[V_h2 < 0.5] <- 0
direcoes <- c(
"Calmo" = 0,
"Nordeste" = 45,
"Leste" = 90,
"Sudeste" = 135,
"Sul" = 180,
"Sudoeste" = 225,
"Oeste" = 270,
"Noroeste" = 315,
"Norte" = 360
)
theta_dir <- sapply(theta2, function(x) {names(direcoes)[which(direcoes == x)]})
theta_dir## [1] "Oeste" "Sul" "Norte" "Leste" "Sudoeste" "Noroeste" "Sudeste"
## [8] "Nordeste" "Calmo"tabela_ex15_b <- data.frame(u, v, V_h2, theta2, theta_dir)
tabela_ex15_b## u v V_h2 theta2 theta_dir
## 1 1.000000e+01 1.836970e-15 10.000 270 Oeste
## 2 -1.224647e-15 1.000000e+01 10.000 180 Sul
## 3 2.449294e-15 -1.000000e+01 10.000 360 Norte
## 4 -1.000000e+01 -6.123234e-16 10.000 90 Leste
## 5 9.999904e+00 9.999904e+00 14.142 225 Sudoeste
## 6 9.999904e+00 -9.999904e+00 14.142 315 Noroeste
## 7 -9.999904e+00 9.999904e+00 14.142 135 Sudeste
## 8 -9.999904e+00 -9.999904e+00 14.142 45 Nordeste
## 9 0.000000e+00 0.000000e+00 0.000 0 CalmoExercício 16
Para as séries de prec_obs e prec_sim calcule: a. A proporção corretamente prevista \(PC= \frac{wc+dc}{n}\) b. O índice de sucesso crítico \(CSI=\frac{wc}{wc+wi+di}\). Onde \(wc\) e \(dc\) são as previsões corretas de dias úmidos \((prec>0.25 \ mm \ dia^{-1})\) e secos respectivamente, \(wi\) e \(di\) são as previsões incorretas de dias úmidos e secos respectivamente. \(n\) é o n° total de previsões.
prec_obs <- c(
0, 0, 0, 0.5, 1, 6, 9, 0.2, 1, 0, 0, 0.25,
10, 15, 8, 3, 0, 0, 0, 0, 0, 0, 0.25, 0,
0, 0, 1, 5, 0, 20, 0, 0, 0, 0, 1, 1,
0, 2, 12, 1, 0, 0, 0, 0, 0, 0, 5, 5
)
prec_sim <- c(
0, 0.2, 0.1, 0, 0, 3, 1, 1, 1, 1, 0, 3,
0, 10, 4, 1, 0.3, 0.5, 0.5, 0.5, 0.5, 0, 0.25, 0.25,
0.25, 0, 0.5, 3, 0, 5, 0, 0, 0, 0, 0.5, 0,
0.25, 0.2, 0, 0.2, 0, 0, 0, 0, 1, 2, 1, 0
)obs_class <- ifelse(prec_obs > 0.25, 1, 0)
sim_class <- ifelse(prec_sim > 0.25, 1, 0)
# Calculando wc, dc, wi, di
wc <- sum(obs_class == 1 & sim_class == 1) # Previsões corretas de dias úmidos
dc <- sum(obs_class == 0 & sim_class == 0) # Previsões corretas de dias secos
wi <- sum(obs_class == 0 & sim_class == 1) # Previsões incorretas de dias úmidos
di <- sum(obs_class == 1 & sim_class == 0) # Previsões incorretas de dias secos
n <- length(prec_obs)
#a)
PC <- (wc + dc) / n
#b)
CSI <- wc / (wc + wi + di)Exercício 17
Escreva o código para calcular as estatísticas abaixo, entre os vetores de valores observados (obs) e previstos (prev) por um dado modelo atmosférico, em um dado local.
- O Viés relativo (%)1. \[PBIAS=100\frac{\sum_{i=1}^{n} (Prev_{i}-Obs_{i})}{\sum_{i=1}^{n} Obs_{i}}\]
- Coeficente de eficiência de Nash-Sutcliffe (NSE)2. \[NSE=1 - \frac{\sum_{i=1}^{n} (Obs_{i}-Sim_{i})²}{\sum_{i=1}^{n} (Obs_{i}-\bar{Obs_{i}})²}\]
v_obs <- c(
-0.49, 0.27, -0.48, 0.8, -1, 0.1, -1.16,
0.58, -1.6, -0.31, 0.45, -0.98, 0.19, 0.73,
-0.49, -0.04, -0.11, 0.46, 2.02, -1.05
)
v_prev <- c(
NA, -0.49, 0.27, -0.48, 0.8, -1, 0.1, -1.16,
0.58, -1.6, -0.31, 0.45, -0.98, 0.19, 0.73,
-0.49, -0.04, -0.11, 0.46, 2.02
)v_obs <- v_obs[!is.na(v_prev)]
v_prev <- v_prev[!is.na(v_prev)]
# a)
PBIAS <- 100 * sum(v_prev - v_obs) / sum(v_obs)
#b)
mean_obs <- mean(v_obs)
NSE <- 1 - sum((v_obs - v_prev)^2) / sum((v_obs - mean_obs)^2)
list(PBIAS = PBIAS, NSE = NSE)## $PBIAS
## [1] -34.5679
##
## $NSE
## [1] -1.691282Mede a tendência média dos valores previstos (ou simulados) em serem maiores (superestimativa) ou menores (subestiva) que os observados. O valor ótimo é 0, menores valores indicam melhor desempenho. Valores positivos indicam tendência de superestimativa e negativos de subestimativa.↩︎
NSE é uma estatística normalizada que que determina a magnitude relativa da variância residual (ruído) comparada a variância dos dados medidos (informação). NSE varia de -Inf a 1. Essencialmente, quanto mais próximo a 1, melhor o modelo.↩︎