Gauss-Seidel Iterative method is essentially the back fitting method for GAM model.
Such as to solve equations
5x1- 2x2+3x3=-1
-3x1+9x2+ x3= 2
2x1 - x2 - 7x3= 3
library(MASS) # package needed for ginv() function
Gauss_Seidel<-function(A,b,x0){
B<-t(A)%*%A
C<-t(A)%*%b
iter<-100 # number of Gauss-Seidel iterations to run,100 iterations
L<-lower.tri(B)*B# lower triang. B,lower.tri() Returns a matrix of logicals,therefore need to *B
U<-upper.tri(B)*B # upper triang. B
D<-diag(diag(B)) # diag of A
n<-1
while(n<=iter){
print(x0)
# Gauss-Seidel formula
x1<-(ginv((L+D)))%*%((-U%*%x0)+C)
x0<-x1
n=n+1
}
return(x0)
}
A<-matrix(c(5,-3,2,-2,9,-1,3,1,-7),ncol=3,nrow=3)
b<-matrix(c(-1,2,3),ncol=1,nrow=3)
x0<-t(t(c(0,0,0)))
Gauss_Seidel(A,b,x0)
After 10 iterations:
> Gauss_Seidel(A,b,x0)
[,1]
[1,] 0
[2,] 0
[3,] 0
[,1]
[1,] -0.1315789
[2,] 0.1380049
[3,] -0.4007323
[,1]
[1,] -0.01103352
[2,] 0.23926763
[3,] -0.41380921
[,1]
[1,] 0.09220576
[2,] 0.28760601
[3,] -0.41850252
[,1]
[1,] 0.1415692
[2,] 0.3105375
[3,] -0.4207159
[,1]
[1,] 0.1649876
[2,] 0.3214148
[3,] -0.4217656
[,1]
[1,] 0.1760960
[2,] 0.3265744
[3,] -0.4222636
[,1]
[1,] 0.1813651
[2,] 0.3290218
[3,] -0.4224998
[,1]
[1,] 0.1838645
[2,] 0.3301827
[3,] -0.4226118
[,1]
[1,] 0.1850501
[2,] 0.3307334
[3,] -0.4226650
[,1]
[1,] 0.1856124
[2,] 0.3309946
[3,] -0.4226902
#Check with solve function.
solve(A,b)
Every time I think I know what's going on, suddenly there's another layer of complications.