Untitled
unknown
r
3 years ago
5.2 kB
8
Indexable
library(ggplot2)
weather <- read.csv("Data/weather.csv")
(Y.model <- lm(rain ~pressure, data = weather))
ggplot(data = weather, aes(x = pressure, y = rain)) + geom_point(size = 1) +
geom_line(aes(y = Y.model$coefficients[1] + pressure*Y.model$coefficients[2]), linetype = "dashed") +
labs(title = "Precipitation (mm) against Air Pressure (hPa)")
# Exponential decay
#INTE FÄRDIGT, cbind(BLANK) ska va ngt annat, körde Y.data, men existerar ej.
press <- data.frame(weather[4])
Y.pred <-
cbind(press,
fit = predict(Y.model, press),
conf = predict(Y.model, press, interval = "confidence"),
pred = predict(Y.model, press, interval = "prediction"))
head(Y.pred)
# get rid of the extra fits
Y.pred$conf.fit <- Y.pred$pred.fit <- NULL
head(Y.pred)
Y.pred$e <- Y.model$residuals
head(Y.pred)
(max.e <- max(abs(Y.pred$e)))
(Y.elims <- c(-max.e, max.e))
## Plot against yhat####
# Add a horizontal line at y=0,
# and expand the y-axis to include +/- max residual.
ggplot(data = Y.pred, aes(x = fit, y = e)) +
geom_point(size = 1) +
geom_hline(yintercept = 0) +
expand_limits(y = Y.elims) +
xlab("Predicted weight loss (g)") +
ylab("Residual") +
labs(tag = "B") +
labs(title = "Residuals vs predicted values Y-hat") +
theme(text = element_text(size = 18))
ggplot(data = weather, aes(x = pressure, y = rain)) + geom_point(size = 1) +
geom_line(aes(y = Y.pred$fit), linetype = "dashed") +
geom_line(aes(y = Y.pred$conf.lwr), linetype = "dashed") +
geom_line(aes(y = Y.pred$conf.upr), linetype = "dashed") +
geom_line(aes(y = Y.pred$pred.lwr), linetype = "dashed") +
geom_line(aes(y = Y.pred$pred.upr), linetype = "dashed") +
labs(title = "Precipitation (mm) against Air Pressure (hPa)")
(Ylog.model <- lm(log(rain) ~ pressure, data = weather))
ggplot(data = weather, aes(x = pressure, y = log(rain))) + geom_point(size = 1) +
geom_line(aes(y = Ylog.model$coefficients[1] + pressure*Ylog.model$coefficients[2]), linetype = "dashed") +
labs(title = "Precipitation (log(mm)) against Air Pressure (hPa)")
# Linear to the left, but after 1010 hPa the data seems to have another slope
## resid anal ###
Ylog.pred <-
cbind(press,
fit = predict(Ylog.model, press),
conf = predict(Ylog.model, press, interval = "confidence"),
pred = predict(Ylog.model, press, interval = "prediction"))
head(Y.pred)
# get rid of the extra fits
Ylog.pred$conf.fit <- Ylog.pred$pred.fit <- NULL
head(Ylog.pred)
Ylog.pred$e <- Ylog.model$residuals
head(Ylog.pred)
(max.e <- max(abs(Ylog.pred$e)))
(Ylog.elims <- c(-max.e, max.e))
## Plot against yhat####
# Add a horizontal line at y=0,
# and expand the y-axis to include +/- max residual.
ggplot(data = Ylog.pred, aes(x = fit, y = e)) +
geom_point(size = 1) +
geom_hline(yintercept = 0) +
expand_limits(y = Ylog.elims) +
xlab("Predicted weight loss (g)") +
ylab("Residual") +
labs(tag = "B") +
labs(title = "Residuals vs predicted values Y-hat") +
theme(text = element_text(size = 18))
ggplot(data = weather, aes(x = pressure, y = rain)) + geom_point(size = 1) +
geom_line(aes(y = exp(Ylog.pred$fit)), linetype = "dashed") +
geom_line(aes(y = exp(Ylog.pred$conf.lwr)), linetype = "dashed") +
geom_line(aes(y = exp(Ylog.pred$conf.upr)), linetype = "dashed") +
geom_line(aes(y = exp(Ylog.pred$pred.lwr)), linetype = "dashed") +
geom_line(aes(y = exp(Ylog.pred$pred.upr)), linetype = "dashed") +
labs(title = "Precipitation (mm) against Air Pressure (hPa)")
(Ycr.model <- lm(rain^(1/3) ~ pressure, data = weather))
ggplot(data = weather, aes(x = pressure, y = rain^(1/3))) + geom_point(size = 1) +
geom_line(aes(y = Ycr.model$coefficients[1] +pressure*Ycr.model$coefficients[2]), linetype = "dashed") +
labs(title = "Precipitation ((mm)^1/3) Against Air Pressure (hPa)")
# Linear, fits our model Y = xB + e
## resid anal ###
Ycr.pred <-
cbind(press,
fit = predict(Ycr.model, press),
conf = predict(Ycr.model, press, interval = "confidence"),
pred = predict(Ycr.model, press, interval = "prediction"))
head(Y.pred)
# get rid of the extra fits
Ycr.pred$conf.fit <- Ycr.pred$pred.fit <- NULL
head(Ycr.pred)
Ycr.pred$e <- Ycr.model$residuals
head(Ycr.pred)
(max.e <- max(abs(Ycr.pred$e)))
(Ycr.elims <- c(-max.e, max.e))
## Plot against yhat####
# Add a horizontal line at y=0,
# and expand the y-axis to include +/- max residual.
ggplot(data = Ycr.pred, aes(x = fit, y = e)) +
geom_point(size = 1) +
geom_hline(yintercept = 0) +
expand_limits(y = Ycr.elims) +
xlab("Predicted weight loss (g)") +
ylab("Residual") +
labs(tag = "B") +
labs(title = "Residuals vs predicted values Y-hat") +
theme(text = element_text(size = 18))
ggplot(data = weather, aes(x = pressure, y = rain)) + geom_point(size = 1) +
geom_line(aes(y = (Ycr.pred$fit)^3), linetype = "dashed") +
geom_line(aes(y = (Ycr.pred$conf.lwr)^3), linetype = "dashed") +
geom_line(aes(y = (Ycr.pred$conf.upr)^3), linetype = "dashed") +
geom_line(aes(y = (Ycr.pred$pred.lwr)^3), linetype = "dashed") +
geom_line(aes(y = (Ycr.pred$pred.upr)^3), linetype = "dashed") +
labs(title = "Precipitation (mm) against Air Pressure (hPa)")
Editor is loading...