Chapter 2: Discrete-Time Markov Chains
Haiyue
19min
Chapter 2: Discrete-Time Markov Chains
2.1 Mathematical Model of Discrete-Time Markov Chains
2.1.1 Rigorous Definition
Let be a stochastic process defined on a discrete time set with state space . If for any and any state sequence , we have:
then is called a discrete-time Markov chain.
2.1.2 Homogeneity Assumption
If the transition probability is independent of time, i.e.:
then it is called a time-homogeneous Markov chain, and the transition probability is called the one-step transition probability.
2.1.3 Complete Characterization of Markov Chains
A homogeneous Markov chain is completely determined by the following two elements:
- Initial Distribution: , where
- Transition Matrix:
2.2 Properties and Computation of Transition Matrices
2.2.1 Basic Properties of Transition Matrix
For a finite state space , the transition matrix is an matrix:
p_{11} & p_{12} & \cdots & p_{1N} \\ p_{21} & p_{22} & \cdots & p_{2N} \\ \vdots & \vdots & \ddots & \vdots \\ p_{N1} & p_{N2} & \cdots & p_{NN} \end{pmatrix}$$ **Properties**: 1. **Non-negativity**: $p_{ij} \geq 0$, $\forall i, j \in S$ 2. **Row-stochasticity**: $\sum_{j=1}^N p_{ij} = 1$, $\forall i \in S$ ### 2.2.2 n-step Transition Probability **Definition**: The probability of going from state $i$ to state $j$ in $n$ steps: $$p_{ij}^{(n)} = P(X_{m+n} = j | X_m = i)$$ **Chapman-Kolmogorov Equation**: $$p_{ij}^{(m+n)} = \sum_{k \in S} p_{ik}^{(m)} p_{kj}^{(n)}$$ **Matrix Form**: $$P^{(n)} = P^n$$ That is, the $n$-step transition matrix equals the $n$-th power of the one-step transition matrix. ### 2.2.3 Evolution of State Distribution Let the state distribution at step $n$ be $\pi^{(n)} = (\pi_1^{(n)}, \pi_2^{(n)}, \ldots)$, where: $$\pi_j^{(n)} = P(X_n = j)$$ Then the evolution of the state distribution follows: $$\pi^{(n+1)} = \pi^{(n)} P$$ Or: $$\pi^{(n)} = \pi^{(0)} P^n$$ ## 2.3 Initial Distribution and Stationary Distribution ### 2.3.1 Impact of Initial Distribution ```mermaid graph LR A[Initial Distribution π⁰] --> B[After one step π¹ = π⁰P] B --> C[After two steps π² = π⁰P²] C --> D[After n steps πⁿ = π⁰Pⁿ] D --> E[Stationary Distribution π = πP] ``` ### 2.3.2 Stationary Distribution **Definition**: If there exists a probability distribution $\pi = (\pi_1, \pi_2, \ldots)$ satisfying: $$\pi = \pi P$$ That is: $$\pi_j = \sum_{i \in S} \pi_i p_{ij}, \quad \forall j \in S$$ then $\pi$ is called the **stationary distribution** of the Markov chain. **Physical Meaning**: If the Markov chain starts from the stationary distribution, then its distribution remains unchanged at any time. ### 2.3.3 Computation of Stationary Distribution The stationary distribution $\pi$ is the solution to the linear system: $$\begin{cases} \pi = \pi P \\ \sum_{i=1}^N \pi_i = 1 \\ \pi_i \geq 0, \quad i = 1, 2, \ldots, N \end{cases}$$ Equivalent to solving: $$\begin{cases} \pi (P - I) = 0 \\ \sum_{i=1}^N \pi_i = 1 \end{cases}$$ ## 2.4 First Passage Time and First Passage Probability ### 2.4.1 First Passage Time **Definition**: The time to first reach state $j$ from state $i$: $$T_{ij} = \inf\{n \geq 1 : X_n = j | X_0 = i\}$$ ### 2.4.2 First Passage Probability **Definition**: The probability of eventually reaching state $j$ from state $i$: $$f_{ij} = P(T_{ij} < \infty | X_0 = i)$$ **Computation Formula**: $$f_{ij} = \sum_{n=1}^{\infty} f_{ij}^{(n)}$$ where $f_{ij}^{(n)}$ is the probability of first reaching state $j$ from state $i$ in exactly $n$ steps. ### 2.4.3 Return Probability and Period - **Return Probability**: $f_{ii}$ represents the probability of eventually returning to state $i$ from state $i$ - **Recurrent State**: If $f_{ii} = 1$, then state $i$ is recurrent - **Transient State**: If $f_{ii} < 1$, then state $i$ is transient ## 2.5 Python Implementation and Applications ```python import numpy as np import matplotlib.pyplot as plt from scipy.linalg import eig import pandas as pd class DiscreteMarkovChain: """Discrete-time Markov Chain Class""" def __init__(self, transition_matrix, state_names=None): """ Initialize discrete-time Markov chain Parameters: ----------- transition_matrix : array-like Transition matrix state_names : list, optional List of state names """ self.P = np.array(transition_matrix, dtype=float) self.n_states = self.P.shape[0] if state_names is None: self.state_names = [f"State_{i}" for i in range(self.n_states)] else: self.state_names = state_names self._validate_transition_matrix() def _validate_transition_matrix(self): """Validate transition matrix""" if self.P.shape[0] != self.P.shape[1]: raise ValueError("Transition matrix must be square") if np.any(self.P < 0): raise ValueError("Transition matrix elements must be non-negative") row_sums = np.sum(self.P, axis=1) if not np.allclose(row_sums, 1.0, rtol=1e-10): raise ValueError("Each row of transition matrix must sum to 1") def n_step_transition_matrix(self, n): """Calculate n-step transition matrix""" if n == 0: return np.eye(self.n_states) elif n == 1: return self.P.copy() else: return np.linalg.matrix_power(self.P, n) def stationary_distribution(self): """Calculate stationary distribution""" # Solve for left eigenvector with eigenvalue 1 eigenvalues, left_eigenvectors = eig(self.P.T) # Find eigenvector corresponding to eigenvalue closest to 1 stationary_idx = np.argmin(np.abs(eigenvalues - 1)) stationary = np.real(left_eigenvectors[:, stationary_idx]) # Normalize to make it a probability distribution stationary = stationary / np.sum(stationary) # Ensure all elements are non-negative if np.any(stationary < 0): stationary = np.abs(stationary) stationary = stationary / np.sum(stationary) return stationary def fundamental_matrix(self): """Calculate fundamental matrix (I - Q)^(-1), where Q is the transient part""" # This assumes all states are transient; in practice, state types should be identified first I = np.eye(self.n_states) try: N = np.linalg.inv(I - self.P) return N except np.linalg.LinAlgError: print("Matrix is not invertible, may have recurrent states") return None def absorption_probabilities(self, transient_states, absorbing_states): """Calculate absorption probabilities""" # Extract Q matrix (transient to transient) and R matrix (transient to absorbing) Q = self.P[np.ix_(transient_states, transient_states)] R = self.P[np.ix_(transient_states, absorbing_states)] # Calculate fundamental matrix I = np.eye(len(transient_states)) N = np.linalg.inv(I - Q) # Calculate absorption probability matrix B = N @ R return B, N def simulate_path(self, initial_state, n_steps, seed=None): """Simulate Markov chain path""" if seed is not None: np.random.seed(seed) if isinstance(initial_state, str): current_state = self.state_names.index(initial_state) else: current_state = initial_state path = [current_state] state_names_path = [self.state_names[current_state]] for _ in range(n_steps): # Choose next state based on current state's transition probabilities next_state = np.random.choice( self.n_states, p=self.P[current_state] ) path.append(next_state) state_names_path.append(self.state_names[next_state]) current_state = next_state return path, state_names_path def expected_return_time(self, state): """Calculate expected return time for a state""" stationary = self.stationary_distribution() if stationary[state] > 0: return 1.0 / stationary[state] else: return np.inf def classify_states(self): """State classification (simplified version)""" # Calculate reachability matrix reachability = np.eye(self.n_states) for i in range(self.n_states): power = np.eye(self.n_states) for _ in range(self.n_states): power = power @ self.P reachability = reachability + power reachability = (reachability > 0).astype(int) # Find strongly connected components (simplified) strongly_connected = [] for i in range(self.n_states): for j in range(self.n_states): if reachability[i, j] and reachability[j, i]: strongly_connected.append((i, j)) return reachability, strongly_connected def plot_transition_graph(self): """Plot state transition graph""" try: import networkx as nx import matplotlib.pyplot as plt G = nx.DiGraph() # Add nodes for i, name in enumerate(self.state_names): G.add_node(i, label=name) # Add edges (only show transitions with probability above threshold) threshold = 0.01 for i in range(self.n_states): for j in range(self.n_states): if self.P[i, j] > threshold: G.add_edge(i, j, weight=self.P[i, j]) # Plot pos = nx.spring_layout(G) plt.figure(figsize=(10, 8)) # Draw nodes nx.draw_networkx_nodes(G, pos, node_color='lightblue', node_size=1000) # Draw edges nx.draw_networkx_edges(G, pos, edge_color='gray', arrows=True, arrowsize=20) # Add labels labels = {i: name for i, name in enumerate(self.state_names)} nx.draw_networkx_labels(G, pos, labels) # Add edge weight labels edge_labels = {(i, j): f"{self.P[i, j]:.2f}" for i in range(self.n_states) for j in range(self.n_states) if self.P[i, j] > threshold} nx.draw_networkx_edge_labels(G, pos, edge_labels) plt.title("Markov Chain State Transition Graph") plt.axis('off') plt.tight_layout() plt.show() except ImportError: print("Need to install networkx library to plot transition graph") # Application Example 1: Random Walk Model def random_walk_example(): """One-dimensional random walk example""" # State space: -2, -1, 0, 1, 2 # Reflection at boundary states n_states = 5 states = list(range(-2, 3)) # [-2, -1, 0, 1, 2] state_names = [f"Position_{i}" for i in states] # Construct transition matrix P = np.zeros((n_states, n_states)) for i in range(n_states): current_pos = states[i] if current_pos == -2: # Left boundary P[i, i] = 0.5 # Stay P[i, i + 1] = 0.5 # Move right elif current_pos == 2: # Right boundary P[i, i - 1] = 0.5 # Move left P[i, i] = 0.5 # Stay else: # Interior states P[i, i - 1] = 0.3 # Move left P[i, i] = 0.4 # Stay P[i, i + 1] = 0.3 # Move right # Create Markov chain rw_chain = DiscreteMarkovChain(P, state_names) # Analyze properties print("Random Walk Markov Chain Analysis") print("=" * 40) # Stationary distribution stationary = rw_chain.stationary_distribution() print("\nStationary distribution:") for i, (name, prob) in enumerate(zip(state_names, stationary)): print(f"{name}: {prob:.4f}") # Expected return time print("\nExpected return time:") for i, name in enumerate(state_names): return_time = rw_chain.expected_return_time(i) print(f"{name}: {return_time:.2f}") # Simulate path path, path_names = rw_chain.simulate_path(2, 20, seed=42) # Start from position 0 print(f"\nSimulated path (20 steps):") print(" -> ".join(path_names[:10]) + " -> ...") return rw_chain # Application Example 2: Brand Loyalty Model def brand_loyalty_example(): """Brand loyalty Markov chain model""" # States: Brand A, Brand B, Brand C state_names = ["Brand_A", "Brand_B", "Brand_C"] # Transition matrix: rows represent current brand, columns represent next purchase brand P = np.array([ [0.7, 0.2, 0.1], # Brand A customers [0.1, 0.8, 0.1], # Brand B customers [0.2, 0.3, 0.5] # Brand C customers ]) brand_chain = DiscreteMarkovChain(P, state_names) print("Brand Loyalty Analysis") print("=" * 30) # Stationary distribution (long-term market share) stationary = brand_chain.stationary_distribution() print("\nLong-term market share:") for name, share in zip(state_names, stationary): print(f"{name}: {share:.1%}") # n-step transition probability print("\n5-step transition probability matrix:") P5 = brand_chain.n_step_transition_matrix(5) df = pd.DataFrame(P5, index=state_names, columns=state_names) print(df.round(4)) # Customer lifecycle analysis print("\nExpected return time (customer loyalty metric):") for i, name in enumerate(state_names): return_time = brand_chain.expected_return_time(i) print(f"{name} customers return to this brand after an average of {return_time:.1f} purchases") return brand_chain # Application Example 3: Income State Transition Model def income_mobility_example(): """Income mobility model""" state_names = ["Low Income", "Middle Income", "High Income"] # Transition matrix based on empirical data P = np.array([ [0.6, 0.3, 0.1], # Low income [0.2, 0.6, 0.2], # Middle income [0.1, 0.4, 0.5] # High income ]) income_chain = DiscreteMarkovChain(P, state_names) print("Income Mobility Analysis") print("=" * 25) # Long-term income distribution stationary = income_chain.stationary_distribution() print("\nLong-term income distribution:") for name, prob in zip(state_names, stationary): print(f"{name}: {prob:.1%}") # Intergenerational mobility analysis print("\nIntergenerational mobility (after 10 years):") P10 = income_chain.n_step_transition_matrix(10) print("Starting state -> State distribution after 10 years") for i, start_state in enumerate(state_names): print(f"{start_state}:") for j, end_state in enumerate(state_names): print(f" -> {end_state}: {P10[i, j]:.1%}") return income_chain # Run all examples if __name__ == "__main__": print("Chapter 2: Discrete-Time Markov Chain Application Examples\n") # Example 1: Random walk rw_chain = random_walk_example() print("\n" + "="*60 + "\n") # Example 2: Brand loyalty brand_chain = brand_loyalty_example() print("\n" + "="*60 + "\n") # Example 3: Income mobility income_chain = income_mobility_example() # Visualization example (requires matplotlib) try: import matplotlib.pyplot as plt # Plot transition probabilities over time fig, axes = plt.subplots(1, 3, figsize=(15, 4)) for idx, (chain, title) in enumerate([ (rw_chain, "Random Walk"), (brand_chain, "Brand Loyalty"), (income_chain, "Income Mobility") ]): steps = range(1, 21) probs = [] for n in steps: Pn = chain.n_step_transition_matrix(n) # Calculate distribution starting from first state prob_dist = Pn[0, :] probs.append(prob_dist) probs = np.array(probs) for i in range(chain.n_states): axes[idx].plot(steps, probs[:, i], label=chain.state_names[i], marker='o') axes[idx].set_title(f"{title}\nState Probability Evolution") axes[idx].set_xlabel("Time Step") axes[idx].set_ylabel("Probability") axes[idx].legend() axes[idx].grid(True, alpha=0.3) plt.tight_layout() plt.show() except ImportError: print("\nNeed to install matplotlib to display graphics") ``` ## 2.6 Convergence Analysis ### 2.6.1 Periodicity **Definition**: The period of state $i$ is: $$d(i) = \gcd\{n \geq 1 : p_{ii}^{(n)} > 0\}$$ - If $d(i) = 1$, state $i$ is called **aperiodic** - If $d(i) > 1$, state $i$ is called **periodic** ### 2.6.2 Limit Theorem For a finite, irreducible, aperiodic Markov chain: $$\lim_{n \to \infty} p_{ij}^{(n)} = \pi_j$$ where $\pi_j$ is the unique stationary distribution. ### 2.6.3 Convergence Rate The convergence rate is determined by the second largest eigenvalue (by modulus) of the transition matrix. If eigenvalues are ordered by modulus as: $$1 = |\lambda_1| > |\lambda_2| \geq |\lambda_3| \geq \cdots$$ then the convergence rate is approximately $|\lambda_2|^n$. ## 2.7 Chapter Summary This chapter provided an in-depth study of the theory and applications of discrete-time Markov chains: **Core Theory**: 1. **Mathematical Model**: Complete characterization of Markov chains through transition matrices 2. **Chapman-Kolmogorov Equation**: Connecting transition probabilities at different time steps 3. **Stationary Distribution**: Probability distribution in long-term steady state 4. **State Classification**: Concepts such as recurrence and periodicity **Key Formulas**: - n-step transition probability: $P^{(n)} = P^n$ - State distribution evolution: $\pi^{(n)} = \pi^{(0)} P^n$ - Stationary distribution: $\pi = \pi P$ - Convergence theorem: $\lim_{n \to \infty} p_{ij}^{(n)} = \pi_j$ **Practical Applications**: 1. **Random Walk**: Spatial diffusion and position models 2. **Brand Loyalty**: Market share prediction 3. **Income Mobility**: Socioeconomic analysis These theories and methods provide an important foundation for financial modeling, particularly in modeling asset price movements, credit state transitions, and other areas with wide applications.