Solving a system of Linear Equations with constraints. Using Symbolic math.

11 views (last 30 days)
I have the following system of equations with a constraint:
n equations with n unknowns
Here are all known. (Its a known n*n matrix of values).
The following is my try at a solution. How can I correct it to get the 's
p_ij = [0.5, 0.5; 0.2, 0.8];
length_pi_i = size(p_ij);
length_pi_i = length_pi_i(1);
syms pi [1 length(p_ij)]
assume(sum(pi) == 1);
zero_vector = zeros(length(pi), 1);
eqn = (pi' - (pi * p_ij)' == zero_vector);
M = solve(eqn, pi);
disp(M);

Answers (2)

Tanmay Das
Tanmay Das on 30 Dec 2021
Hi,
I tried reproducing your steps in MATLAB R2021b release and received the following output in the command window:
pi1: 2/7
pi2: 5/7
I tried verifying the answer manually and the answer seems correct to me.
In case you want to learn more about symbolic functions, you may look into 'syms','assume', and 'solve' function documentation.

John D'Errico
John D'Errico on 2 Jan 2023
Edited: John D'Errico on 2 Jan 2023
This is just an eigenvalue problem. (Yes, I know that is probably beyond the scope of the question that was asked. But this is just a basic Markov chain problem. And we should recognize the matrix P_ij as a 2x2 Markov transition matrix.)
p_ij = [0.5, 0.5; 0.2, 0.8]
p_ij = 2×2
0.5000 0.5000 0.2000 0.8000
But first, we can ask if a solution of the form requested exists for this problem? That is, does the problem
X*P_ij = X
Where X is a vector of length 2, have a solution? For that to be the case, then we must have 1 as an eigenvalue of this matrix.
format long g
[V,D] = eig(p_ij.')
V = 2×2
-0.707106781186547 -0.371390676354104 0.707106781186547 -0.928476690885259
D = 2×2
0.3 0 0 1
So in fact, it is true that 1 is an eigenvalue. We see that because we see that D(2,2) == 1. The corresponding vector, thus V(:,2) is the solution.
pi_vec = V(:,2)'
pi_vec = 1×2
-0.371390676354104 -0.928476690885259
But we should also know that an eigenvector can be scaled arbitrarily. We can multiply it by any constant, and it is still an equally valid eigenvector. So we can now just rescale pi_vec to sum to 1 as needed. And by switching to format rat, we find the desired long term probabilities for each of the two states of this Markov process.
format rat
pi_vec = pi_vec/sum(pi_vec)
pi_vec =
2/7 5/7
There was not even a need to use solve. Just eig was sufficient here. Could it have been done using syms? Of course, though I refuse to name a variable pi, as that tramples on the terribly useful number pi.
syms X [1,2]
T = sym(p_ij) %% the transition matrix, in symbolic form.
T = 
solve(X*T == X,sum(X) == 1,X)
ans = struct with fields:
X1: 2/7 X2: 5/7
Again, we see the expected solution drops out.

Products


Release

R2021a

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!