There are 3 ways using matrices and linear algebra :
Method number 1 - Gauss-Jordan Elimination
Write the linear system in the form Ax=bAx=b where
A=[(4,-1),(5,-3)] is the coefficient matrix of the system,
x=[(x),(y)] is the column vector of variables, and
b=[(-3),(-23)] is the column vector of solutions.
Now augment the coefficient matrix with the column vector of solutions and perform elementary row operations until you obtain the 2xx2 identity matrix I_2 on the left hand side and then what remains on the right hand side will be the solution set for x. I will attach a sketch at the end to show these steps.
Eventually we get :
[(4,-1),(5,-3)]|[(-3),(-23)] -> [(1,0),(0,1)]|[(2),(11)]
and hence x=2 and y=11.
Method number 2 - Using the inverse matrix method
Since Ax=b it implies that x=A^(-1)b.
We may find A^(-1) by using elementary row operations in augmenting the coefficient matrix A with the identity matrix I_2, alternatively, we may use the formula :
A^(-1)=1/(det(A))*adj(A), where det(A)=Delta.
=1/(-12-(-5))*[(-3,1),(-5,4)]
=[(3/7,-1/7),(5/7,-4/7)]
Hence x=A^(-1)b=[(3/7,-1/7),(5/7,-4/7)][(-3),(-23)]=[(2),(11)]
and so x=2 and y=11.
Method number 3 - Kramer's Rule
Find the determinants of the matrices formed by replacing each column vector in the coefficient matrix with the column vector of solutions :
Delta_x=|(-3,-1),(-23,-3)|=(-3)(-3)-(-1)(-23)=-14.
Delta_y=|(4,-3),(5,-23)|=-77.
Then : x=Delta_x/Delta=-14/-7=2
and y=Delta_y/Delta=-77/-7=11.