HW1 S102B
Isaiah C. Mireles
Link to website :
https://rpubs.com/Isaiah-Mireles/1420507
Q1)
a)
cos_theta <-
function(v1,v2){
as.numeric(v1%*%v2)/(sqrt(sum((v1)^2))*sqrt(sum((v2)^2)))
}
# problem with direction :
# get_angle <- function(v1,v2){acos(cos_theta(v1, v2))*(180/pi)}
get_angle <- function(v1, v2){
angle1 <- atan2(v1[2], v1[1]) * (180/pi) # put each component of vec
angle2 <- atan2(v2[2], v2[1]) * (180/pi)
# compare angle of vec
angle_diff <- angle2 - angle1
# ensure its 0 to 360
angle_diff <- (angle_diff + 360) %% 360
angle_diff
}b)
Okay, so this question is about the rotation matrix
\[ R(\theta) =\begin{bmatrix}\cos\theta & -\sin\theta \\\sin\theta & \cos\theta\end{bmatrix} \]
unit_vec <- c(1,0) # pointing right
rotate_vector <- function(numeric_vector, degrees){
# check dimension
if( length(numeric_vector) == 2 ){
# convert degrees to radians
theta <- degrees * (pi/180)
# create rotation matrix
rotation_matrix <- matrix(c(cos(theta), -sin(theta),
sin(theta), cos(theta)),
nrow = 2,
byrow = T)
# rotate vector
rotated_vector <- rotation_matrix %*% numeric_vector |> as.vector()
# check angle using get_angle : how much did it rotate?
# copare to i-unit vector pointing right ccw
angle_check <- get_angle(c(1,0), rotated_vector)
return(list(rotated_vector = rotated_vector,
angle_check = angle_check))
} else {
print("Vector must be 2D")
}
}example <-
rotate_vector(unit_vec, 90)$rotated_vector
rotate_vector(unit_vec, 90)## $rotated_vector
## [1] 6.123234e-17 1.000000e+00
##
## $angle_check
## [1] 9017 degrees
rotate_vector(unit_vec, 17)## $rotated_vector
## [1] 0.9563048 0.2923717
##
## $angle_check
## [1] 17110 degrees
rotate_vector(unit_vec, 110)## $rotated_vector
## [1] -0.3420201 0.9396926
##
## $angle_check
## [1] 110270 degrees
rotate_vector(unit_vec, 270) ## $rotated_vector
## [1] -1.83697e-16 -1.00000e+00
##
## $angle_check
## [1] 270Q2)
Q3)
Attempt 1)
okay so A is square and so is P
P is orth
\[ \text{Def. } P_{\text{Orth}} \implies P'P=P'P=I \\ \implies P^{-1}=P' \]
all of the columns are orthogonal^
Each feature is a unit vector^
Preserves length & angles under rotation^
an eigen value and vector are
\[ \text{Def. } \lambda \vec{v}=A\vec{v} \]
so in other words, there exists a matrix which scales and/or flips the vector.
eigen vectors are the directions which stay the same after transformation in a matrix
- an eigen value is the amt of scaling
–
Attempt 2)
\[ \text{Let } A \in \mathbb{R}^{k \times k} \text{ and let } P \in \mathbb{R}^{k \times k} \text{ be orthogonal, so } P^T = P^{-1}. \]
\[ \text{Define } M = P^T A P. \]
\[ \text{Since } P^T = P^{-1}, \text{ we can write} \quad M = P^{-1} A P, \]
\[ \text{so } M \text{ is similar to } A! \]
\[ \text{Let } \lambda \text{ be an eigenvalue of } A \text{ with eigenvector } x \neq 0, \\ \text{ so} \quad Ax = \lambda x. \]
\[ \text{Define } y = P^T x. \]
\[ \text{Since } P \text{ is invertible, } y \neq 0. \]
\[ \text{Now compute:} \quad My = (P^T A P) y = P^T A P (P^T x). \]
\[ \text{Using } PP^T = I, \quad My = P^T A x. \]
\[ \text{Substitute } Ax = \lambda x: \quad My = P^T (\lambda x) = \lambda P^T x = \lambda y. \]
\[ \text{Thus, } y \text{ is an eigenvector of } M \text{ with eigenvalue } \lambda. \]
\[ \text{Therefore, every eigenvalue of } A \text{ is an eigenvalue of } P^T A P. \]
\[ \text{By symmetry (reversing the argument), they have the same eigenvalues.} \]
Q4)
a)
- its actually \(-\sqrt{3}\) for off diagnols
r1 <- c(5, -sqrt(3))
r2 <- c(-sqrt(3), 7)
m <- matrix(c(r1,r2), nrow=2, byrow = T)
m## [,1] [,2]
## [1,] 5.000000 -1.732051
## [2,] -1.732051 7.000000det <- det(m); det## [1] 32trace <- sum(diag(m)); trace## [1] 12eigen_val <- (trace + sqrt(trace^2 - 4*det(m))*c(1,-1))/2
eigen_val## [1] 8 4m - 4*diag(1,2)## [,1] [,2]
## [1,] 1.000000 -1.732051
## [2,] -1.732051 3.000000
Graphic)
b)
so as we can see that our eigen vectors are just gonna be a scalar multiple of what we got before :
by adjusting \(\frac{c}{\sqrt{\lambda}}\) we obtain the
correctly scaled vectors :
Solution :
\[ \text{For : } c^2=1 \\ \implies \\ \boxed{ \vec{v}_1=<\frac{\sqrt{3}}{4},\frac{1}{4}> \\ \text{and,} \\ \vec{v}_2 = <\frac{1}{4\sqrt{2}},-\frac{1}{4\sqrt{2}}> } \]
\[ \text{For : } c^2=\pi \\ \implies \\ \boxed{ \vec{v}_1= <\frac{\pi\sqrt{3}}{4},\frac{\pi}{4}>\\ \text{and,} \\ \vec{v}_2 = <\frac{\pi}{4\sqrt{2}},-\frac{\pi \sqrt{3}}{4\sqrt{2}}> } \]
- notice we just mult the \(c^2=1\) example by \(\pi\) , showing its just a scalar multiple as stated previously.
Q6)
a)
Recall :
\[ \hat{\beta} = (X'X)^{-1}X'\vec{y} \]
Where \(X\) is our design matrix : \(X=[\vec{1} \ \vec{x}]\)
# getwd()
df <- read.csv("data1.csv")Therefore :
intercept <- rep(1, nrow(df))
design_mat <- cbind(intercept, df$x)
y <- df$y
beta_hat <-
solve(t(design_mat)%*%design_mat)%*%t(design_mat)%*%y
mdl <- lm(y~x, data=df)
summary(mdl)$coefficients[,1]## (Intercept) x
## -4.269400 4.485172beta_hat## [,1]
## intercept -4.269400
## 4.485172- as we see, we get the same result ^
b)
plot(df$x, df$y)
abline(mdl, col = "red")
- least squares linear mdl doesn’t describe behavior
c)
x <- df$x
# ?poly
model <- lm(y ~ poly(x, degree = 6, raw = T), data = df)
model2 <- lm(y ~ x + I(x^2) + I(x^3) + I(x^4) + I(x^5) + I(x^5) + I(x^6), data = df)- personal note : These output the same thing – however, using
poly()function – the preset israw = Fso it makes each regressor orthogonal. This is not the same as just the other version usingI()
print("Polyn Reg")## [1] "Polyn Reg"summary(model)$r.squared## [1] 0.8254print("Linear Reg")## [1] "Linear Reg"summary(mdl)$r.squared## [1] 0.4777358- it appears to be doing much better, but i bet its overfit.
d)
plot(x,y)
lines(smooth.spline(x, predict(model)), col = "red")
abline(mdl, col = "blue")
- it definitely describes the behavior of the data better but i worry that with a 6th order polynomial we may be fitting on to noise, especially on the tail ends
library(tidymodels)
library(recipes)
# 1. Create LOOCV
loocv_folds <- vfold_cv(df, v = nrow(df)) # folds = sample sz
# 2. Define recipe with a tunable degree for a polynomial
poly_rec <- recipe(y ~ x, data = df) %>%
step_poly(x, degree = tune()) # Tune the degree
# 3. Model specification
poly_spec <- linear_reg() %>% set_engine("lm")
# 4. Workflow
poly_wf <- workflow() %>%
add_recipe(poly_rec) %>%
add_model(poly_spec)
# 5. Tune grid (testing degrees 1 through 6)
degree_grid <- grid_regular(degree(range = c(1, 6)), levels = 6)
tune_results <- tune_grid(
poly_wf,
resamples = loocv_folds,
grid = degree_grid
)
# 6. Select best model
best_poly <- select_best(tune_results, metric = "rmse")mdl_metrics <-
collect_metrics(tune_results) |> filter(.metric=="rmse") |> select(degree, .metric, mean, std_err)
mdl_metrics## # A tibble: 6 × 4
## degree .metric mean std_err
## <dbl> <chr> <dbl> <dbl>
## 1 1 rmse 12.3 0.973
## 2 2 rmse 12.3 0.841
## 3 3 rmse 8.51 0.455
## 4 4 rmse 8.60 0.468
## 5 5 rmse 8.60 0.435
## 6 6 rmse 8.62 0.439- it appears the best performing model is degree 3
poly3_mdl <- lm(y ~ poly(x, degree = 3, raw = T), data = df)
poly6_mdl <- lm(y ~ poly(x, degree = 6, raw = T), data = df)
plot(x,y)
lines(smooth.spline(x, predict(poly3_mdl)), col = "red")
lines(smooth.spline(x, predict(poly6_mdl)), col = "green")
abline(mdl, col = "blue")
Linear regression : BLUE
Polynomial regression – degree 3 : RED
Polynomial regression – degree 6 : GREEN
mdl_metrics[c(3,6),]## # A tibble: 2 × 4
## degree .metric mean std_err
## <dbl> <chr> <dbl> <dbl>
## 1 3 rmse 8.51 0.455
## 2 6 rmse 8.62 0.439- here we notice a dramatic change in
std_errandmean(rsme) when degree 3 hits
paste0("poly3_mdl :",summary(poly3_mdl)$r.squared)## [1] "poly3_mdl :0.820499460621341"paste0("poly6_mdl :",summary(poly6_mdl)$r.squared)## [1] "poly6_mdl :0.825399979221013"- further notice despite a dramatically higher degree we only see about half a percent better in \(r^2\)
round(summary(poly6_mdl)$r.squared - summary(poly3_mdl)$r.squared , 3)## [1] 0.005Personal Graphing function :
incomplete!!!!
need to fix radians, degrees output
does handle above 180 degrees properly
graph_vector2D <- function(v1, v2) {
# Determine plot limits
x_limits <- range(c(0, v1[1], v2[1]))
y_limits <- range(c(0, v1[2], v2[2]))
# Create an empty plot
plot(
NA,
xlim = x_limits,
ylim = y_limits,
xlab = "X-axis",
ylab = "Y-axis",
main = "2D Vector Plot",
asp = 1
)
# Add grid and axes
grid()
abline(h = 0, v = 0, col = "gray")
# Draw vectors from origin
arrows(
x0 = 0, y0 = 0,
x1 = v1[1], y1 = v1[2],
col = "blue", lwd = 2, length = 0.1
)
arrows(
x0 = 0, y0 = 0,
x1 = v2[1], y1 = v2[2],
col = "red", lwd = 2, length = 0.1
)
# Add labels near the vector tips
text(v1[1], v1[2], labels = "v1", pos = 4, col = "blue")
text(v2[1], v2[2], labels = "v2", pos = 4, col = "red")
}
# use case :
graph_vector2D(-c(0,1),c(-1,0))
graph_vector3D <- function(v1, v2) {
# Load required package
if (!requireNamespace("plotly", quietly = TRUE)) {
stop("The 'plotly' package is required. Please install it using install.packages('plotly').")
}
# Helper function to add arrowheads using cones
create_arrowhead <- function(plot, vec, color) {
plotly::add_trace(
plot,
type = "cone",
x = vec[1],
y = vec[2],
z = vec[3],
u = vec[1],
v = vec[2],
w = vec[3],
sizemode = "absolute",
sizeref = 0.4,
anchor = "tip",
showscale = FALSE,
colorscale = list(c(0, color), c(1, color)),
hoverinfo = "none"
)
}
# Create base plot with vector shafts
plot <- plotly::plot_ly() |>
plotly::add_trace(
x = c(0, v1[1]),
y = c(0, v1[2]),
z = c(0, v1[3]),
type = "scatter3d",
mode = "lines+text",
line = list(color = "blue", width = 6),
text = c("", "v1"),
textposition = "top center",
name = "v1"
) |>
plotly::add_trace(
x = c(0, v2[1]),
y = c(0, v2[2]),
z = c(0, v2[3]),
type = "scatter3d",
mode = "lines+text",
line = list(color = "red", width = 6),
text = c("", "v2"),
textposition = "top center",
name = "v2"
)
# Add arrowheads correctly
plot <- create_arrowhead(plot, v1, "blue")
plot <- create_arrowhead(plot, v2, "red")
# Add layout settings
plot <- plotly::layout(
plot,
title = "3D Vector Plot",
scene = list(
xaxis = list(title = "X-axis"),
yaxis = list(title = "Y-axis"),
zaxis = list(title = "Z-axis")
)
)
return(plot)
}
# use case :
# Plot using your function
graph_vector3D(c(1, 0, 0), c(0, 1, 0))graph_vector <- function(v1, v2,
table_position = "topright",
docu = FALSE) {
# -------------------------
# Documentation Mode
# -------------------------
if (docu) {
cat("
graph_vector(v1, v2, table_position = 'topright', docu = FALSE)
Description:
Plots two vectors in 2D or 3D (dimension must be < 4).
Automatically computes and displays:
- cos(theta)
- theta in radians
- theta in degrees
Parameters:
v1 Numeric vector (length 2 or 3).
v2 Numeric vector (same dimension as v1).
table_position Character string specifying where the
angle information appears.
Possible values:
'topright'
'topleft'
'bottomleft'
'bottomright'
docu Logical (TRUE/FALSE).
If TRUE, prints this documentation only.
Output:
- 2D vectors → Static base R plot.
- 3D vectors → Interactive Plotly plot.
Examples:
graph_vector(c(1,0), c(0,1))
graph_vector(c(1,0,0), c(0,1,0), table_position='bottomleft')
")
return(invisible(NULL))
}
# -------------------------
# Validation
# -------------------------
if (!is.numeric(v1) || !is.numeric(v2)) {
stop("Both v1 and v2 must be numeric vectors.")
}
if (length(v1) != length(v2)) {
stop("v1 and v2 must have the same dimension.")
}
dim <- length(v1)
if (dim >= 4) {
stop("Vectors must have dimension less than 4 (2D or 3D only).")
}
valid_positions <- c(
"topright", "topleft",
"bottomleft", "bottomright"
)
if (!(table_position %in% valid_positions)) {
stop("Invalid table_position value.")
}
# -------------------------
# Compute angle info
# -------------------------
cos_val <- cos_theta(v1, v2)
angle_rad <- get_angle(v1, v2)
angle_deg <- angle_rad * 180 / pi
info_text <- c(
paste0("cos(theta) = ", round(cos_val, 4)),
paste0("theta (rad) = ", round(angle_rad, 4)),
paste0("theta (deg) = ", round(angle_deg, 2))
)
# -------------------------
# 2D Case
# -------------------------
if (dim == 2) {
graph_vector2D(v1, v2)
legend(
table_position,
legend = info_text,
bty = "n"
)
# -------------------------
# 3D Case
# -------------------------
} else if (dim == 3) {
plot <- graph_vector3D(v1, v2)
pos_map <- list(
topright = list(x = 1, y = 1),
topleft = list(x = 0, y = 1),
bottomleft = list(x = 0, y = 0),
bottomright = list(x = 1, y = 0)
)
coords <- pos_map[[table_position]]
plot <- plotly::layout(
plot,
annotations = list(
list(
x = coords$x,
y = coords$y,
xref = "paper",
yref = "paper",
text = paste(info_text, collapse = "<br>"),
showarrow = FALSE,
align = "left"
)
)
)
return(plot)
}
}# Use cases :
graph_vector(c(0, -1), c(-1, 0))
graph_vector(
c(2, 3, 1),
c(3, -2, 0),
table_position = "bottomleft"
)v1 <- c(1, 0)
v2 <- c(0, 1)
graph_vector(v1, v2, table_position = "topright")
Documentation example :
graph_vector(docu = TRUE)##
## graph_vector(v1, v2, table_position = 'topright', docu = FALSE)
##
## Description:
## Plots two vectors in 2D or 3D (dimension must be < 4).
## Automatically computes and displays:
## - cos(theta)
## - theta in radians
## - theta in degrees
##
## Parameters:
## v1 Numeric vector (length 2 or 3).
## v2 Numeric vector (same dimension as v1).
## table_position Character string specifying where the
## angle information appears.
##
## Possible values:
## 'topright'
## 'topleft'
## 'bottomleft'
## 'bottomright'
##
## docu Logical (TRUE/FALSE).
## If TRUE, prints this documentation only.
##
## Output:
## - 2D vectors → Static base R plot.
## - 3D vectors → Interactive Plotly plot.
##
## Examples:
## graph_vector(c(1,0), c(0,1))
## graph_vector(c(1,0,0), c(0,1,0), table_position='bottomleft')Negative cosine examples :
# Base vector
y1 <- c(1, 0)
# Opposite direction
y2 <- c(-1, 0)
# Mostly opposite
y3 <- c(-1, 1)
# Another opposite-ish
y4 <- c(-2, -1)
# Slightly opposite
y5 <- c(-0.5, 0.2)
graph_vector(y1, y2)
graph_vector(y1, y3)
graph_vector(y1, y4, table_position="bottomright")
graph_vector(y1, y5)
Visualize Q2B)
graph_vector(unit_vec, example)
graph_vector(unit_vec, rotate_vector(unit_vec, 270)$rotated_vector)