<SYSTEM_TASK:>
Given the following text description, write Python code to implement the functionality described below step by step
<END_TASK>
<USER_TASK:>
Description:
 State Value Function
Next we will find different ways to compute the Value function given by a stochastic policy $\pi(s) = p(a\mid s)$.
We want to calculate $V_{\pi}(s)$ (the state-value-function given a policy)
Here we draw an Markov Decision Process (MDP) with three states $\mathcal{S}={s_1,s_2,s_3}$ and three possible actions $\mathcal{A}={a_1,a_2,a_3}$, moving to state $s_1$, moving to state $s_2$ and moving to state $s_3$.
Step1:  Policy Evaluation by Dynamic Programming
For the MDP represented above we define the state transition probability matrix $\mathcal{P}^a_{ss'}=p(S_{t+1}=s'\mid S_{t}=s, A_t=a)$. In this MDP we assume that when we choose to move to state $s_i$, $i={1,2,3}$ we always end up in that state, meaning that $\mathcal{P}^a_{ss'}=p(S_{t+1}=s'\mid S_{t}=s, A_t=a)=1$. In this case, $\mathcal{P}^{\pi}=\mathcal{P}^a_{ss'}\pi(a\mid s) = \pi(a\mid s)$ the Bellman Expectation equation becomes (Check page 14 and 16 from the lecture slides.)
Step2:  Policy Evaluation by Linear Programming
The state-value-function can be directly solved through linear programming (as shown on page 15 from the lecture slides)
Step3:  The result stays the same. 
Policy Evaluation by Monte Carlo Sampling
We can design yet another way of evaluating the value of a given policy $\pi$, see lecture slides pag.20.
The intuition is to incrementally the expected return from sampled episodes, sequences of triplets ${(s_i,a_i,r_{i})}{i=1}^N$. The function $\color{blue}{gt}$ computes the total discounted reward from a list of sequential rewards obtained by sampling the policy
Step4:  As can be seen the result is nearly the same as the state-value-function calculated above.
So far we have seen different ways of given a known policy $\pi(a\mid s)$ how to comput its value $V_{\pi}(s)$. Next, we wish to find the optimal policy $\pi^\ast(s)$ for the MDP in the example.
Policy Optimization by Q-Learning
This code solves a very easy problem
Step5:  Value Iteration
<ASSISTANT_TASK:>
Python Code:
%load_ext autoreload
%autoreload 2
Explanation: State Value Function
Next we will find different ways to compute the Value function given by a stochastic policy $\pi(s) = p(a\mid s)$.
We want to calculate $V_{\pi}(s)$ (the state-value-function given a policy)
Here we draw an Markov Decision Process (MDP) with three states $\mathcal{S}={s_1,s_2,s_3}$ and three possible actions $\mathcal{A}={a_1,a_2,a_3}$, moving to state $s_1$, moving to state $s_2$ and moving to state $s_3$.
End of explanation
import numpy as np
policy = np.array([[0.3, 0.2, 0.5], [0.5, 0.4, 0.1], [0.8, 0.1, 0.1]])
print("This is represents the policy with 3 states and 3 actions p(row=a|col=s):\n", np.matrix(policy))
# 'raw_rewards' variable contains rewards obtained after transition to each state
# In our example it doesn't depend on source state
raw_rewards = np.array([1.5, -1.833333333, 19.833333333])
# 'rewards' variable contains expected values of the next reward for each state
rewards = np.matmul(policy, raw_rewards)
assert np.allclose(rewards, np.array([10., 2., 3.]))
gamma = 0.1
print('This are the rewards for each action:\n', rewards)
state_value_function = np.array([0 for i in range(3)])
print('Policy evaluation:')
for i in range(20):
    print('V_{}={}'.format(i, state_value_function))
    state_value_function = rewards + gamma * (np.matmul(policy, state_value_function))
print('\nV={}'.format(state_value_function))
Explanation: Policy Evaluation by Dynamic Programming
For the MDP represented above we define the state transition probability matrix $\mathcal{P}^a_{ss'}=p(S_{t+1}=s'\mid S_{t}=s, A_t=a)$. In this MDP we assume that when we choose to move to state $s_i$, $i={1,2,3}$ we always end up in that state, meaning that $\mathcal{P}^a_{ss'}=p(S_{t+1}=s'\mid S_{t}=s, A_t=a)=1$. In this case, $\mathcal{P}^{\pi}=\mathcal{P}^a_{ss'}\pi(a\mid s) = \pi(a\mid s)$ the Bellman Expectation equation becomes (Check page 14 and 16 from the lecture slides.):
$$
V_{\pi}(s) = \sum_{a\in\mathcal{A}} \pi(a\mid s)\left( \mathcal{R}^a_s + \gamma \sum_{s'\in \mathcal{S}}\mathcal{P}^a_{ss'}V_{\pi}(s')\right) = \mathcal{R}^{\pi}+ \gamma \sum_{s'\in \mathcal{S}}\pi(a\mid s)V_{\pi}(s')
$$
End of explanation
solution=np.matmul(np.linalg.inv(np.eye(3)-0.1*policy), rewards)
print('Solution by inversion:\nV={}'.format(state_value_function))
Explanation: Policy Evaluation by Linear Programming
The state-value-function can be directly solved through linear programming (as shown on page 15 from the lecture slides):
$$
V_{\pi}(s)=\left(I-\gamma\mathcal{P}^{\pi}\right)^{-1}\mathcal{R}^{\pi}
$$
End of explanation
import random
from collections import defaultdict
reward_counter = np.array([0., 0., 0.])
visit_counter = np.array([0., 0., 0.])
nIterations = 400
def gt(rewardlist, gamma=0.1):
    '''
    Function to calculate the total discounted reward
    >>> gt([10, 2, 3], gamma=0.1)
    10.23
    '''
    total_disc_return = 0
    for (i, value) in enumerate(rewardlist):
        total_disc_return += (gamma ** i) * value
    return total_disc_return
for i in range(nIterations): 
    start_state = random.randint(0, 2)
    next_state = start_state
    rewardlist = []
    occurence = defaultdict(list)
    for i in range(250): #draw samples from the policy recursively over horizon of N=250
        rewardlist.append(rewards[next_state])
        occurence[next_state].append(len(rewardlist) - 1)
        action = np.random.choice(np.arange(0, 3), p=policy[next_state])
        next_state = action
    for state in occurence:
        for value in occurence[state]: #update state value function E[G_t|s]=S(s)/N(s)
            rew = gt(rewardlist[value:]) 
            reward_counter[state] += rew # S(s)
            visit_counter[state] += 1 # N(s)
print("MC policy evaluation V=", reward_counter / visit_counter)
Explanation: The result stays the same. 
Policy Evaluation by Monte Carlo Sampling
We can design yet another way of evaluating the value of a given policy $\pi$, see lecture slides pag.20.
The intuition is to incrementally the expected return from sampled episodes, sequences of triplets ${(s_i,a_i,r_{i})}{i=1}^N$. The function $\color{blue}{gt}$ computes the total discounted reward from a list of sequential rewards obtained by sampling the policy: $G_t=r_t+\gamma r{t+1}+\gamma^2 r_{t+2}+\dots+\gamma^N r_{t+N}$.
The value of a policy can also be computed by looking at its empirical expected cumulative discounted return:
$$
V_{\pi}(s) = \mathbb{E}_{\pi}\left[G_t\mid S_t=s\right]
$$
End of explanation
q_table = np.zeros((3, 3)) #state action value function Q-table
gamma = 0.1
alpha = 1.0
eps = 0.1
def get_eps_greedy_action(state):
    if random.uniform(0, 1) < eps:
        return random.randint(0, 2)
    return np.argmax(q_table[state]).item()
for i in range(1001):
    state = random.randint(0, 2)
    action = get_eps_greedy_action(state)
    next_state = action
    reward = raw_rewards[next_state]
    next_q = max(q_table[next_state]) #s.a. value evaluation at the next state
    q_table[state, action] = q_table[state, action] + alpha* (
            reward + gamma * (next_q) - q_table[state, action])  #Q-Table update
    if i % 200 == 0:
        print("Q_{}(s,a)=".format(i),q_table)
Explanation: As can be seen the result is nearly the same as the state-value-function calculated above.
So far we have seen different ways of given a known policy $\pi(a\mid s)$ how to comput its value $V_{\pi}(s)$. Next, we wish to find the optimal policy $\pi^\ast(s)$ for the MDP in the example.
Policy Optimization by Q-Learning
This code solves a very easy problem: using the rewards it calculates the optimal action-value-function (page 26 on slides).
It samples a state-action pair randomly, so that all state-action pairs can be seen.
End of explanation
import numpy as np
raw_rewards = np.array([1.5, -1.833333333, 19.833333333])
gamma = 0.1
state_value_function = np.zeros(3)
print('V_{} = {}'.format(0, state_value_function))
for i in range(1000):
    for s in range(3):
        Q_s = [raw_rewards[s_next] + gamma * state_value_function[s_next]
               for s_next in range(3)]
        state_value_function[s] = max(Q_s)
    if i % 100 == 99:
        print('V_{} = {}'.format(i + 1, state_value_function))
Explanation: Value Iteration
End of explanation
<END_TASK>