LeastSquaresDegree3 <- function(t, b) { # Usage: t is independent variable vector, b is function data # i.e., b = f(t) # This code only works for finding second degree polynomial best-fit models t <- t(t) dataSize <- length(b) A <- mat.or.vec(dataSize, 4) # Built-in R function to create zero matrix or zero vector of arbitrary size # Given basis phi(z) = 1 + z + z^2 + z^3 # Define matrix A A[1:dataSize,1] = 1 A[1:dataSize,2] = t A[1:dataSize,3] = t^2 A[1:dataSize,4] = t^3 # Compute QR decomposition of A. Pull Q and R out of QRdecomp QRdecomp <- qr(A) Q <- qr.Q(QRdecomp, complete=TRUE) R <- qr.R(QRdecomp, complete=TRUE) # Perform Q^* b^T (Conjugate transpose of Q) c <- Conj(t(Q))%*%t(b) # Find x. R isn't square - so we have to use qr.solve x <- qr.solve(R, c) cat("coefficients are:", x[1],",",x[2],",",x[3],",",x[4],".") # return(c(x=x)) # Now plot it. Least squares "l" plot first, then the points in red. plot(t, x[1] + x[2]*t + x[3]*t^2 + x[4]*t^3, type='l', xlab="independent variable t", ylab="function values f(t)", main="Data Plotted with 3rd Degree Least Squares Polynomial", col="blue") points(t, b, col="red") # Code by Eric Johnson -- Focus | Numeric }