# PageRank

#### Librerías utilizadas


```python
import numpy as np
import math
import pandas as pd
from scipy.sparse import csr_matrix 
  
# Para las pruebas
from pytest import approx
import pprint

# Para los grafos
import networkx as nx
import matplotlib.pyplot as plt
```

## Definimos una set para practicar (por ahora 4 páginas)

Donde podemos observar en la columnas:
- `pagina`: la página a la que hace referencia.
- `outlinks`: los links a los que cada `pagina` hace referencia.


```python
pages_links = pd.DataFrame(columns=['pagina','outlinks'])
```


```python
pages_links['pagina'] = ['A', 'B', 'C', 'D']
```


```python
pages_links['outlinks'][0] = ['C', 'B']
```


```python
pages_links['outlinks'][1] = ['D']
```


```python
pages_links['outlinks'][2] = ['A', 'B', 'D']
```


```python
pages_links['outlinks'][3] = ['C']
```


```python
pages_links
```




<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }

    .dataframe tbody tr th {
        vertical-align: top;
    }

    .dataframe thead th {
        text-align: right;
    }
</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th>pagina</th>
      <th>outlinks</th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>0</th>
      <td>A</td>
      <td>[C, B]</td>
    </tr>
    <tr>
      <th>1</th>
      <td>B</td>
      <td>[D]</td>
    </tr>
    <tr>
      <th>2</th>
      <td>C</td>
      <td>[A, B, D]</td>
    </tr>
    <tr>
      <th>3</th>
      <td>D</td>
      <td>[C]</td>
    </tr>
  </tbody>
</table>
</div>



## Definimos un grafo pequeño de acuerdo a los datos anteriores


```python
nodos = pages_links['pagina']

plt.figure(figsize=(8,8))
G = nx.DiGraph()

for nodo in nodos:
    G.add_node(nodo)
    
for page, outlinks in zip(pages_links['pagina'],pages_links['outlinks']):
    for outlink in outlinks:
        G.add_edge(page, outlink)

pos = nx.spring_layout(G)


nx.draw_networkx_nodes(G, pos, node_color='#a0a88a', node_size=400)
nx.draw_networkx_labels(G, pos)
nx.draw_networkx_edges(G, pos, edge_color='#305029', width=1.8, arrows = True, style='dashed')

plt.show()
```

### Definimos nuestra matriz estocástica conforme a lo visto en el grafo:


```python
sparseMatrix = csr_matrix((len(pages_links), len(pages_links)),  
                          dtype = np.float).toarray() 
```


```python
newdf = pd.DataFrame(sparseMatrix, columns = pages_links['pagina'], index=pages_links['pagina'])
```


```python
count = 0
for page, outlinks in zip(pages_links['pagina'],pages_links['outlinks']):
    score = 1/len(outlinks)
    dict_page = {page: np.zeros(len(pages_links))}
    for outlink in outlinks:
        newdf[page][outlink] = score
```


```python
newdf
```




<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }

    .dataframe tbody tr th {
        vertical-align: top;
    }

    .dataframe thead th {
        text-align: right;
    }
</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th>pagina</th>
      <th>A</th>
      <th>B</th>
      <th>C</th>
      <th>D</th>
    </tr>
    <tr>
      <th>pagina</th>
      <th></th>
      <th></th>
      <th></th>
      <th></th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>A</th>
      <td>0.0</td>
      <td>0.0</td>
      <td>0.333333</td>
      <td>0.0</td>
    </tr>
    <tr>
      <th>B</th>
      <td>0.5</td>
      <td>0.0</td>
      <td>0.333333</td>
      <td>0.0</td>
    </tr>
    <tr>
      <th>C</th>
      <td>0.5</td>
      <td>0.0</td>
      <td>0.000000</td>
      <td>1.0</td>
    </tr>
    <tr>
      <th>D</th>
      <td>0.0</td>
      <td>1.0</td>
      <td>0.333333</td>
      <td>0.0</td>
    </tr>
  </tbody>
</table>
</div>




```python
H = newdf.to_numpy()
```


```python
pd.DataFrame(newdf)
```




<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }

    .dataframe tbody tr th {
        vertical-align: top;
    }

    .dataframe thead th {
        text-align: right;
    }
</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th>pagina</th>
      <th>A</th>
      <th>B</th>
      <th>C</th>
      <th>D</th>
    </tr>
    <tr>
      <th>pagina</th>
      <th></th>
      <th></th>
      <th></th>
      <th></th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>A</th>
      <td>0.0</td>
      <td>0.0</td>
      <td>0.333333</td>
      <td>0.0</td>
    </tr>
    <tr>
      <th>B</th>
      <td>0.5</td>
      <td>0.0</td>
      <td>0.333333</td>
      <td>0.0</td>
    </tr>
    <tr>
      <th>C</th>
      <td>0.5</td>
      <td>0.0</td>
      <td>0.000000</td>
      <td>1.0</td>
    </tr>
    <tr>
      <th>D</th>
      <td>0.0</td>
      <td>1.0</td>
      <td>0.333333</td>
      <td>0.0</td>
    </tr>
  </tbody>
</table>
</div>



## Probamos la obtención de Eigenvectores con Factorización QR (algoritmo desarrollado en la práctica 1.2)


```python
def sign(x):
    """
    Helper function for computing sign of real number x.
    """
    if x >= 0:
        return 1
    else:
        return -1
```


```python
def inverse_matrix(A):
    """
    Helper function for computing the inverse of matrix A nxn.
    """
    inv_temp = []
    
    for i in range (0,A.shape[0]): 
        e1 = np.zeros((A.shape[0],1))
        e1[i] = 1
        A_inv_i = np.linalg.solve(A.T, e1)
        inv_temp.append(A_inv_i)  
        
    for j in range (0,len(inv_temp)-1):
        if j == 0:
            r_inv = np.row_stack((inv_temp[j].T, inv_temp[j+1].T))
        else:
            r_inv = np.row_stack((r_inv, inv_temp[j+1].T))
        
    return(r_inv)
```


```python
def off(T):
    """
    Helper function for computing off(T).
    """
    mask = np.ones(T.shape, dtype = bool)
    np.fill_diagonal(mask, 0)
    T_sin_dig = T[mask]
    sum = 0
    for i in range (0,len(T_sin_dig)):
        # print (T_sin_dig[i])
        sum += T_sin_dig[i]**2
    off_T = np.sqrt(sum)
    
    return(off_T)
```


```python
def QR_Householder(A):
    """
    Decompose a real square matrix A where: A = QR usign Householder reflections.
    Args:
        A (numpy ndarray): Matrix in which QR algorithm will be performed.
    Returns:
        Q (numpy ndarray): orthogonal matrix (its columns are orthogonal unit vectors meaning
                            Q.TQ = QQ.T = I)
        R (numpy ndarray): upper triangular matrix.
    """
    a_aux = []
    
    for j in range (0,A.shape[0]-1):          #Asumimos que A SIEMPRE será cuadrada
        e = np.zeros((A.shape[0]-j, 1))
        e[0] = 1
        if j == 0:
            v = A[:,0] + sign(A[0,0])*np.linalg.norm(A[:,0])*e.T[0]
            beta = 2/v.dot(v)
            aux = A[:,0:] - beta*np.outer(v,v.dot(A[:,0:]))
            a_aux.append(aux)
        else:
            v = aux[1:,1] + sign(aux[1,1])*np.linalg.norm(aux[1:,1])*e.T[0]
            beta = 2/v.dot(v)  
            aux = aux[1:,1:] - beta*np.outer(v,v.dot(aux[1:,1:]))
            a_aux.append(aux)
    
    n,n = A.shape
    count = 0
    
    if n == 2:
        R = a_aux[0]
    else:
        for h in range(n-2,0,-1):
            if count == 0:
                temp_col = np.column_stack((np.zeros(2 + count), a_aux[h]))
            else:
                temp_col = np.column_stack((np.zeros(2 + count), temp_row))
            temp_row = np.row_stack((a_aux[h-1][0, 0:], temp_col))
            count += 1
        R = temp_row
    R_inv = inverse_matrix(R)
    
    Q = A@R_inv   
    
    return (Q,R)
```


```python
def Eig_Householder(A):
    """
    Compute QR algorithm usign Householder reflections multiple times to approximate eigenvalues and 
    eigenvectors of A.
    Args:
        A (numpy ndarray): Symmetric matrix 
    Returns:
        evalues (numpy ndarray): Array with eigenvalues of matrix A.
        evectors (numpy ndarray): Matrix with eigenvectors of matrix A on its columns.
    """
    tol = 10e-8
    max_iters = 2000
    tk_fro_norm = np.linalg.norm(A, 'fro')
    iterations = 0
    off_Tk = off(A)
    
    q_aux = []
    
    while off_Tk > tol*tk_fro_norm and iterations < max_iters:
        iterations += 1
        if iterations == 1:
            Tk_minus1 = A
        else:
            Tk_minus1 = Tk
        
        q,r = QR_Householder(Tk_minus1)    
        q_aux.append(q)
        
        Tk = r@q
        
        off_Tk = off(Tk)
        tk_fro_norm = np.linalg.norm(Tk, 'fro')
    
    #Multiplicación Q0*Q1...Qn para cálculo de eigenvectores
    for g in range(0,len(q_aux)-1,1):
        if g == 0:
            q_temp = q_aux[g]@q_aux[g+1]
        else:
            q_temp = q_temp@q_aux[g+1]
    
    if iterations == 0:
        evalues = "error"
        evectors = "error"
    else:
        evectors = q_temp
        evalues = np.diag(Tk)
    
    return (evalues,evectors)
```


```python
eigval, eigvec = Eig_Householder(H)
```


```python
pd.DataFrame(eigval)
```




<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }

    .dataframe tbody tr th {
        vertical-align: top;
    }

    .dataframe thead th {
        text-align: right;
    }
</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th>0</th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>0</th>
      <td>1.000000</td>
    </tr>
    <tr>
      <th>1</th>
      <td>-0.626538</td>
    </tr>
    <tr>
      <th>2</th>
      <td>-0.357205</td>
    </tr>
    <tr>
      <th>3</th>
      <td>-0.016257</td>
    </tr>
  </tbody>
</table>
</div>




```python
pd.DataFrame(eigvec)
```




<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }

    .dataframe tbody tr th {
        vertical-align: top;
    }

    .dataframe thead th {
        text-align: right;
    }
</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th>0</th>
      <th>1</th>
      <th>2</th>
      <th>3</th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>0</th>
      <td>0.232495</td>
      <td>0.529635</td>
      <td>0.237150</td>
      <td>-0.780508</td>
    </tr>
    <tr>
      <th>1</th>
      <td>0.348743</td>
      <td>0.192872</td>
      <td>0.785518</td>
      <td>0.473434</td>
    </tr>
    <tr>
      <th>2</th>
      <td>0.697486</td>
      <td>-0.672574</td>
      <td>0.004505</td>
      <td>-0.247260</td>
    </tr>
    <tr>
      <th>3</th>
      <td>0.581238</td>
      <td>0.479512</td>
      <td>-0.571577</td>
      <td>0.324855</td>
    </tr>
  </tbody>
</table>
</div>



#### No olvidemos comprobar que los eigenvalores y eigenvectores sean los esperados


```python
eigval_eig , eigvec_eig = np.linalg.eig(H)
```


```python
pd.DataFrame(eigval_eig)
```




<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }

    .dataframe tbody tr th {
        vertical-align: top;
    }

    .dataframe thead th {
        text-align: right;
    }
</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th>0</th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>0</th>
      <td>1.000000+0.000000j</td>
    </tr>
    <tr>
      <th>1</th>
      <td>-0.626538+0.000000j</td>
    </tr>
    <tr>
      <th>2</th>
      <td>-0.186731+0.480774j</td>
    </tr>
    <tr>
      <th>3</th>
      <td>-0.186731-0.480774j</td>
    </tr>
  </tbody>
</table>
</div>




```python
pd.DataFrame(eigvec_eig)
```




<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }

    .dataframe tbody tr th {
        vertical-align: top;
    }

    .dataframe thead th {
        text-align: right;
    }
</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th>0</th>
      <th>1</th>
      <th>2</th>
      <th>3</th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>0</th>
      <td>-0.232495+0.000000j</td>
      <td>0.445747+0.000000j</td>
      <td>0.146497+0.377184j</td>
      <td>0.146497-0.377184j</td>
    </tr>
    <tr>
      <th>1</th>
      <td>-0.348743+0.000000j</td>
      <td>0.090025+0.000000j</td>
      <td>0.435929+0.112414j</td>
      <td>0.435929-0.112414j</td>
    </tr>
    <tr>
      <th>2</th>
      <td>-0.697486+0.000000j</td>
      <td>-0.837833+0.000000j</td>
      <td>-0.626088+0.000000j</td>
      <td>-0.626088-0.000000j</td>
    </tr>
    <tr>
      <th>3</th>
      <td>-0.581238+0.000000j</td>
      <td>0.302061+0.000000j</td>
      <td>0.043661-0.489599j</td>
      <td>0.043661+0.489599j</td>
    </tr>
  </tbody>
</table>
</div>



## Probamos con un Método de la Potencia adaptado a PageRank (obtenido de [Global Software Support](https://www.youtube.com/watch?v=kSmQbVxqOJc&feature=emb_title))

Este método de la potencia se define como: 
    
$$ PR_{t+1} = H PR_t $$

Donde:
- $H$: la matriz estocástica.
- $PR$: PageRank de la página indicada.

#### Definimos un vector inicial 

El vector será de tamaño $n$ (número de nodos) y la inicialización de cada valor será $1/n$.

En este caso:

\begin{align}
    v &= \begin{bmatrix}
           1/4 \\
           1/4 \\
           1/4 \\
           1/4
         \end{bmatrix}
  \end{align}


```python
v = np.ones(len(pages_links))
```


```python
v = v*1/len(v)
```


```python
pprint.pprint(v)
```

    array([0.25, 0.25, 0.25, 0.25])


Para obtener el segundo vector se debe hacer:
    
$$v_2 = Hv$$


```python
v_2 = H@v

pprint.pprint(v_2)
```

    array([0.08333333, 0.20833333, 0.375     , 0.33333333])


Para obtener el tercer vector:
    
$$ v_3 = Hv_2 = H(Hv) = H^2v$$


```python
v_3  = H@v_2

pprint.pprint(v_3)
```

    array([0.125     , 0.16666667, 0.375     , 0.33333333])



```python
v_3_r = (H@H)@v

pprint.pprint(v_3_r)
```

    array([0.125     , 0.16666667, 0.375     , 0.33333333])


Entonces, observamos que se puede representar de la forma: 
$$v_n = H^{n}v$$

Si iteramos varias veces tiende a un valor equilibrado:


```python
n = 20
H = np.linalg.matrix_power(H, n)
v_n = H@v
```


```python
v_n
```




    array([0.12500368, 0.18750089, 0.37499297, 0.31250246])




```python
df = pd.DataFrame(v_n, columns=['pagerank'])
```


```python
df['nodo'] = nodos
```


```python
df = df.set_index('nodo')
```


```python
power_method = df.pagerank.to_dict()
```


```python
power_method
```




    {'A': 0.1250036838777434,
     'B': 0.18750089306127105,
     'C': 0.37499296714248953,
     'D': 0.3125024559184955}




```python
power_method_round = {} 

for node in power_method:
    if node in power_method:
        power_method_round[node] = round(power_method[node],2)
```


```python
power_method_round
```




    {'A': 0.13, 'B': 0.19, 'C': 0.37, 'D': 0.31}



## Probamos con el Método de la Potencia visto en clase


```python
n = 4
q_k = v
max_iter = 20
q_k_iter = np.zeros((n, max_iter))
lambda_k_iter = np.zeros(max_iter)

for k in range(max_iter):
    z_k = H@q_k
    q_k = z_k/np.linalg.norm(z_k)
    q_k_iter[:,k] = q_k
    lambda_k = q_k.T@H@q_k
    lambda_k_iter[k] = lambda_k
```


```python
pd.DataFrame(q_k)
```




<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }

    .dataframe tbody tr th {
        vertical-align: top;
    }

    .dataframe thead th {
        text-align: right;
    }
</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th>0</th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>0</th>
      <td>0.232495</td>
    </tr>
    <tr>
      <th>1</th>
      <td>0.348743</td>
    </tr>
    <tr>
      <th>2</th>
      <td>0.697486</td>
    </tr>
    <tr>
      <th>3</th>
      <td>0.581238</td>
    </tr>
  </tbody>
</table>
</div>



## Ahora probamos con PageRank de `Networkx`

El cálculo del eigenvector se realiza mediante el método iterativo de la potencia y no tiene garantía de convergencia. Las iteraciones se detendrán después de que se alcance una tolerancia de error de `len(G)*tol`. Si el número de iteraciones excede el número máximo de iteraciones permitido, se genera una excepción `networkx.exception.PowerIterationFailedConvergence`.
El algoritmo PageRank fue diseñado para grafos dirigidos, pero este algoritmo no verifica si el grafo de entrada es dirigido o no, y se ejecutará en grafos no dirigidos, convirtiendo cada borde en el grafo dirigido en dos bordes.


```python
from networkx.algorithms.link_analysis.pagerank_alg import pagerank
```


```python
pr = pagerank(G, max_iter=20)
```


```python
pr
```




    {'A': 0.13867229558179425,
     'B': 0.1976082981211456,
     'C': 0.3570799577952133,
     'D': 0.30663944850184677}




```python
pagerank_round = {}    

for node in pr:
    if node in pr:
        pagerank_round[node] = round(pr[node],2)
```


```python
pagerank_round
```




    {'A': 0.14, 'B': 0.2, 'C': 0.36, 'D': 0.31}



## Ahora comparamos todos los resultados obtenidos


```python
pagerank_round_pd = pd.DataFrame.from_dict(pagerank_round, orient='index')
```


```python
eigvector_comparison = pd.DataFrame({"NX PageRank": pagerank_round_pd[0],
                                    "Potencia PR": v_n,
                                    "Potencia OPT-2020": q_k,
                                    "LinAlg": eigvec_eig[:,0].real,
                                    "QR_Eq4": eigvec[:,0]}
                                    ).reset_index().iloc[:,1:]
eigvector_comparison
```




<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }

    .dataframe tbody tr th {
        vertical-align: top;
    }

    .dataframe thead th {
        text-align: right;
    }
</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th>NX PageRank</th>
      <th>Potencia PR</th>
      <th>Potencia OPT-2020</th>
      <th>LinAlg</th>
      <th>QR_Eq4</th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>0</th>
      <td>0.14</td>
      <td>0.125004</td>
      <td>0.232495</td>
      <td>-0.232495</td>
      <td>0.232495</td>
    </tr>
    <tr>
      <th>1</th>
      <td>0.20</td>
      <td>0.187501</td>
      <td>0.348743</td>
      <td>-0.348743</td>
      <td>0.348743</td>
    </tr>
    <tr>
      <th>2</th>
      <td>0.36</td>
      <td>0.374993</td>
      <td>0.697486</td>
      <td>-0.697486</td>
      <td>0.697486</td>
    </tr>
    <tr>
      <th>3</th>
      <td>0.31</td>
      <td>0.312502</td>
      <td>0.581238</td>
      <td>-0.581238</td>
      <td>0.581238</td>
    </tr>
  </tbody>
</table>
</div>



En la tabla anterior podemos apreciar que tanto el eigenvector obtenido por `Networkx` como el obtenido por el método de la potencia adaptado son muy parecidos. De hecho, éstos son factores de los vectores obtenidos con los otros 3 métodos. Adicionalmente, se aprecia que el orden se mantiene, por lo que se confirma que los PageRanks obtenidos son los correctos.

#### Comparamos método de la potencia vs `pagerank_alg.pagerank`


```python
assert pagerank_round == approx(power_method_round, rel=2e-1)
```

## Mismo grafo con el  PageRank obtenido por el método de la potencia


```python
plt.figure(figsize=(8,8))
G = nx.DiGraph()

for nodo in power_method:
    G.add_node(nodo)

for page, outlinks in zip(pages_links['pagina'],pages_links['outlinks']):
    for outlink in outlinks:
        G.add_edge(page, outlink)

pos = nx.spring_layout(G)

labels = {}    
for node in G.nodes():
    if node in power_method:
        #set the node name as the key and the label as its value 
        labels[node] = "{}:{}".format(node,round(power_method[node],2))

# Es lo mismo que power_method pero redondeado        
print(labels)

nx.draw_networkx_nodes(G, pos, node_color='#a0a88a', node_size=400)
nx.draw_networkx_labels(G, pos, labels=labels, font_weight="bold")
nx.draw_networkx_edges(G, pos, edge_color='#305029', width=1.8, arrows = True, style='dashed')

plt.show()
```

Entonces podemos observar que la página con mayor puntuación es la `C`.

## Ahora intentaremos con otro set de datos (37 páginas)

[set de datos de Stanford](http://snap.stanford.edu/data/web-Stanford.html)


```python
path = "../data/web-Stanford.txt"
```


```python
def read_csv_Stanford(path, nrow=103):
    """
    Read the txt given, only works for web-Stanford 
    """
    df = pd.read_csv(path, nrows=nrow)#, usecols=['FromNode', 'ToNode'])
    # Remove the index
    df = df[3:]
    df = df.rename(columns = {'# Directed graph (each unordered pair of nodes is saved once): web-Stanford.txt ': 
                         'FromTo'})
    df[['FromNode','ToNode']] = df.FromTo.str.split('\t',expand=True) 
    df = df.drop('FromTo', axis=1)
    return df
```


```python
df_stanford = read_csv_Stanford(path,103)
```


```python
df_stanford
```




<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }

    .dataframe tbody tr th {
        vertical-align: top;
    }

    .dataframe thead th {
        text-align: right;
    }
</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th>FromNode</th>
      <th>ToNode</th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>3</th>
      <td>1</td>
      <td>6548</td>
    </tr>
    <tr>
      <th>4</th>
      <td>1</td>
      <td>15409</td>
    </tr>
    <tr>
      <th>5</th>
      <td>6548</td>
      <td>57031</td>
    </tr>
    <tr>
      <th>6</th>
      <td>15409</td>
      <td>13102</td>
    </tr>
    <tr>
      <th>7</th>
      <td>2</td>
      <td>17794</td>
    </tr>
    <tr>
      <th>...</th>
      <td>...</td>
      <td>...</td>
    </tr>
    <tr>
      <th>98</th>
      <td>74361</td>
      <td>226411</td>
    </tr>
    <tr>
      <th>99</th>
      <td>74361</td>
      <td>239874</td>
    </tr>
    <tr>
      <th>100</th>
      <td>74361</td>
      <td>255694</td>
    </tr>
    <tr>
      <th>101</th>
      <td>74361</td>
      <td>264216</td>
    </tr>
    <tr>
      <th>102</th>
      <td>74361</td>
      <td>266083</td>
    </tr>
  </tbody>
</table>
<p>100 rows × 2 columns</p>
</div>




```python
def transform_data():
    """
    Transform the nodes that are not in the FromNode and ToNode match
    This only works for 103 nrows
    """
    missing_from = list(set(df_stanford['FromNode']) - set(df_stanford['ToNode']) )
    missing_to = list(set(df_stanford['ToNode'])  - set(df_stanford['FromNode']) )
    mitad = int(len(missing_to)/2)
    for val in missing_to[:mitad]:
        df_stanford['ToNode'][df_stanford['ToNode'] == val] = missing_from[0]
    for val in missing_to[mitad:]:
        df_stanford['ToNode'][df_stanford['ToNode'] == val] = missing_from[1]
```


```python
transform_data()
```


```python
# Hacemos una agrupación y la lista de los outlinks
paginas_stanford = df_stanford.groupby('FromNode').aggregate(lambda tdf: tdf.unique().tolist()).reset_index()
```


```python
paginas_stanford = paginas_stanford.rename(columns={'FromNode': 'pagina',
                                                    'ToNode': 'outlinks'})
```


```python
paginas_stanford.head(10)
```




<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }

    .dataframe tbody tr th {
        vertical-align: top;
    }

    .dataframe thead th {
        text-align: right;
    }
</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th>pagina</th>
      <th>outlinks</th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>0</th>
      <td>1</td>
      <td>[6548, 15409]</td>
    </tr>
    <tr>
      <th>1</th>
      <td>100193</td>
      <td>[2]</td>
    </tr>
    <tr>
      <th>2</th>
      <td>102355</td>
      <td>[2]</td>
    </tr>
    <tr>
      <th>3</th>
      <td>105318</td>
      <td>[2]</td>
    </tr>
    <tr>
      <th>4</th>
      <td>105730</td>
      <td>[2, 3]</td>
    </tr>
    <tr>
      <th>5</th>
      <td>115926</td>
      <td>[2]</td>
    </tr>
    <tr>
      <th>6</th>
      <td>140864</td>
      <td>[2]</td>
    </tr>
    <tr>
      <th>7</th>
      <td>15409</td>
      <td>[1]</td>
    </tr>
    <tr>
      <th>8</th>
      <td>163550</td>
      <td>[2]</td>
    </tr>
    <tr>
      <th>9</th>
      <td>164599</td>
      <td>[2]</td>
    </tr>
  </tbody>
</table>
</div>




```python
paginas_stanford.shape
```




    (37, 2)




```python
def plot_node_graph(dataframe_pages, power_method=None, with_labels=False):
    """
    Plotea el grafo
    :param[in] dataframe_pages: dataframe de nx2
                El dataframe debe contener una columna llamada página y otra llamada outlinks;
                donde outlinks contenga una lista de los outlinks por página
    """
    nodos = dataframe_pages['pagina']

    plt.figure(figsize=(10,10))
    G = nx.DiGraph()

    for nodo in nodos:
        G.add_node(nodo)

    for page, outlinks in zip(dataframe_pages['pagina'], dataframe_pages['outlinks']):
        for outlink in outlinks:
            G.add_edge(page, outlink)

    pos = nx.spring_layout(G)


    nx.draw_networkx_nodes(G, pos, node_color='#a0a88a', node_size=400)
    if with_labels:
        labels = {}    
        for node in G.nodes():
            if node in power_method:
                #set the node name as the key and the label as its value 
                labels[node] = "{}:{}".format(node,round(power_method[node],2))
        nx.draw_networkx_labels(G, pos, labels=labels, font_weight="bold")
    else:
        nx.draw_networkx_labels(G, pos)
    nx.draw_networkx_edges(G, pos, edge_color='#305029', width=1.8, arrows = True, style='dashed')

    plt.show()
    return G
```


```python
G37 = plot_node_graph(paginas_stanford)
```

#### Definimos nuestra matriz estocástica


```python
def get_h_matrix(dataframe_pages):
    """
    Obtiene la matriz H, la matriz estocástica
    :param[in] dataframe_pages: dataframe de nx2
                El dataframe debe contener una columna llamada página y otra llamada outlinks;
                donde outlinks contenga una lista de los outlinks por página
    """
    sparseMatrix = csr_matrix((len(dataframe_pages), len(dataframe_pages)),  
                          dtype = np.float).toarray() 
    newdf = pd.DataFrame(sparseMatrix, columns = dataframe_pages['pagina'], index=dataframe_pages['pagina'])
    count = 0
    for page, outlinks in zip(dataframe_pages['pagina'],dataframe_pages['outlinks']):
        score = 1/len(outlinks)
        dict_page = {page: np.zeros(len(dataframe_pages))}
        for outlink in outlinks:
            newdf[page][outlink] = score
    H = newdf.to_numpy()
    return H
```


```python
H37 = get_h_matrix(paginas_stanford)
```


```python
H37.shape
```




    (37, 37)




```python
from networkx.algorithms.link_analysis.pagerank_alg import google_matrix
```


```python
matrix_37 = google_matrix(G37, nodelist=paginas_stanford.pagina)
```


```python
A_37 = matrix_37.T
```

### Obtenemos PageRank por el Método de la Potencia adaptado a PR


```python
def power_method_pagerank(dataframe_pages, H, iteraciones):
    """
    Obtiene el pagerank 
    :param[in] dataframe_pages: dataframe de nx2
                El dataframe debe contener una columna llamada página y otra llamada outlinks;
                donde outlinks contenga una lista de los outlinks por página
    :param[in] H: Matriz estocástica
    :param[in] iteraciones: número de iteraciones a realizar
    
    :return power_method: Un diccionario con la página como llave y el pagerank como valor
    """
    # Obtenemos nuestro vector inicial
    v = np.ones(len(dataframe_pages))
    v = v*1/len(v)
    
    # Iteramos
    H = np.linalg.matrix_power(H, iteraciones)
    v_n = H@v
    
    # Convertimos a dataframe
    print(v_n.shape)
    df = pd.DataFrame(v_n.T, columns=['pagerank'])
    
    df['nodo'] = dataframe_pages['pagina']
    
    df = df.set_index('nodo')
    
    power_method = df.pagerank.to_dict()
    
    return power_method
```


```python
power_method_dic = power_method_pagerank(paginas_stanford, A_37, iteraciones=100)
```

    (1, 37)



```python
power_method_dic
```




    {'1': 0.10093628803480034,
     '100193': 0.013330403674657317,
     '102355': 0.013330403674657317,
     '105318': 0.013330403674657317,
     '105730': 0.013330403674657317,
     '115926': 0.013330403674657317,
     '140864': 0.013330403674657317,
     '15409': 0.04695200491257369,
     '163550': 0.013330403674657317,
     '164599': 0.013330403674657317,
     '175799': 0.013330403674657317,
     '17794': 0.013330403674657317,
     '178642': 0.013330403674657317,
     '181714': 0.013330403674657317,
     '190453': 0.013330403674657317,
     '2': 0.33831392139907474,
     '204189': 0.013330403674657317,
     '204604': 0.013330403674657317,
     '210870': 0.013330403674657317,
     '213966': 0.013330403674657317,
     '225119': 0.013330403674657317,
     '241596': 0.013330403674657317,
     '243294': 0.013330403674657317,
     '246897': 0.013330403674657317,
     '251658': 0.013330403674657317,
     '25202': 0.013330403674657317,
     '252915': 0.013330403674657317,
     '280935': 0.013330403674657317,
     '3': 0.026783358060477675,
     '53625': 0.013330403674657317,
     '54582': 0.013330403674657317,
     '64930': 0.013330403674657317,
     '6548': 0.04695200491257369,
     '73764': 0.013330403674657317,
     '74361': 0.026819908766137674,
     '84477': 0.013330403674657317,
     '98628': 0.013330403674657317}




```python
# Redondeamos a 2 dígitos
power_method_round = {}    

for node in power_method_dic:
    if node in power_method_dic:
        #set the node name as the key and the label as its value 
        power_method_round[node] = round(power_method_dic[node],2)
```


```python
G37_pot = plot_node_graph(paginas_stanford, power_method=power_method_round, with_labels=True)
```

No se puede apreciar pero:


```python
MaxKey = max(power_method_round, key=power_method_round.get)
```


```python
print("la página con más rank es: {} con un rank de {}".format(MaxKey, power_method_round[MaxKey]))
```

    la página con más rank es: 2 con un rank de 0.34


### Obtenemos PageRank con el paquete `Networkx`


```python
PR_stanford = pagerank(G37, max_iter=1000)
```


```python
PR_round_stanford = {} 

for node in PR_stanford:
    if node in PR_stanford:
        PR_round_stanford[node] = round(PR_stanford[node],2)
```


```python
MaxKeyPR = max(PR_round_stanford, key=PR_round_stanford.get)
```


```python
print("La página con más rank es: {} con un rank de {}".format(MaxKeyPR, PR_round_stanford[MaxKeyPR]))
```

    La página con más rank es: 2 con un rank de 0.34


## Obtenemos con el Método de la Potencia visto en clase


```python
def power_method_opt_pagerank(dataframe_pages, H, iteraciones):
    """
    Obtiene el pagerank 
    :param[in] dataframe_pages: dataframe de nx2
                El dataframe debe contener una columna llamada página y otra llamada outlinks;
                donde outlinks contenga una lista de los outlinks por página
    :param[in] H: Matriz estocástica
    :param[in] iteraciones: número de iteraciones a realizar
    
    :return power_method: Un diccionario con la página como llave y el pagerank como valor
    """
    # Obtenemos nuestro vector inicial
    v = np.ones(len(dataframe_pages))
    v = v*1/len(v)
    
    n = len(v)
    q_k = v
    q_k = q_k.reshape(n,1)
    max_iter = iteraciones

    for k in range(max_iter):
        z_k = H@q_k
        q_k = z_k/np.linalg.norm(z_k)
        lambda_k = q_k.T @ H @ q_k
    
    # Convertimos a dataframe
    df = pd.DataFrame(q_k, columns=['pagerank'])
    
    df['nodo'] = dataframe_pages['pagina']
    
    df = df.set_index('nodo')
    
    power_method_opt = df.pagerank.to_dict()
    
    return power_method_opt
```


```python
power_method_opt_dictionary = power_method_opt_pagerank(paginas_stanford, A_37, iteraciones=100)
```


```python
MaxKey_pot_opt = max(power_method_opt_dictionary, key=power_method_opt_dictionary.get)
```


```python
print("La página con más rank es: {} con un rank de {}".format(MaxKey_pot_opt, power_method_opt_dictionary[MaxKey_pot_opt]))
```

    La página con más rank es: 2 con un rank de 0.9173865700385999


## Comparando resultados


```python
eigvector_comparison_37 = pd.DataFrame({"Página": PR_round_stanford.keys(),
                                    "NX PageRank": PR_round_stanford.values(),
                                    "Potencia PR": power_method_round.values(),
                                    "Potencia OPT-2020": power_method_opt_dictionary.values()}
                                    #"LinAlg": eigvec_eig[:,0].real}
                                    ).reset_index().iloc[:,1:]
eigvector_comparison_37
```




<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }

    .dataframe tbody tr th {
        vertical-align: top;
    }

    .dataframe thead th {
        text-align: right;
    }
</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th>Página</th>
      <th>NX PageRank</th>
      <th>Potencia PR</th>
      <th>Potencia OPT-2020</th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>0</th>
      <td>1</td>
      <td>0.10</td>
      <td>0.10</td>
      <td>0.273703</td>
    </tr>
    <tr>
      <th>1</th>
      <td>100193</td>
      <td>0.01</td>
      <td>0.01</td>
      <td>0.036147</td>
    </tr>
    <tr>
      <th>2</th>
      <td>102355</td>
      <td>0.01</td>
      <td>0.01</td>
      <td>0.036147</td>
    </tr>
    <tr>
      <th>3</th>
      <td>105318</td>
      <td>0.01</td>
      <td>0.01</td>
      <td>0.036147</td>
    </tr>
    <tr>
      <th>4</th>
      <td>105730</td>
      <td>0.01</td>
      <td>0.01</td>
      <td>0.036147</td>
    </tr>
    <tr>
      <th>5</th>
      <td>115926</td>
      <td>0.01</td>
      <td>0.01</td>
      <td>0.036147</td>
    </tr>
    <tr>
      <th>6</th>
      <td>140864</td>
      <td>0.01</td>
      <td>0.01</td>
      <td>0.036147</td>
    </tr>
    <tr>
      <th>7</th>
      <td>15409</td>
      <td>0.05</td>
      <td>0.05</td>
      <td>0.127317</td>
    </tr>
    <tr>
      <th>8</th>
      <td>163550</td>
      <td>0.01</td>
      <td>0.01</td>
      <td>0.036147</td>
    </tr>
    <tr>
      <th>9</th>
      <td>164599</td>
      <td>0.01</td>
      <td>0.01</td>
      <td>0.036147</td>
    </tr>
    <tr>
      <th>10</th>
      <td>175799</td>
      <td>0.01</td>
      <td>0.01</td>
      <td>0.036147</td>
    </tr>
    <tr>
      <th>11</th>
      <td>17794</td>
      <td>0.01</td>
      <td>0.01</td>
      <td>0.036147</td>
    </tr>
    <tr>
      <th>12</th>
      <td>178642</td>
      <td>0.01</td>
      <td>0.01</td>
      <td>0.036147</td>
    </tr>
    <tr>
      <th>13</th>
      <td>181714</td>
      <td>0.01</td>
      <td>0.01</td>
      <td>0.036147</td>
    </tr>
    <tr>
      <th>14</th>
      <td>190453</td>
      <td>0.01</td>
      <td>0.01</td>
      <td>0.036147</td>
    </tr>
    <tr>
      <th>15</th>
      <td>2</td>
      <td>0.34</td>
      <td>0.34</td>
      <td>0.917387</td>
    </tr>
    <tr>
      <th>16</th>
      <td>204189</td>
      <td>0.01</td>
      <td>0.01</td>
      <td>0.036147</td>
    </tr>
    <tr>
      <th>17</th>
      <td>204604</td>
      <td>0.01</td>
      <td>0.01</td>
      <td>0.036147</td>
    </tr>
    <tr>
      <th>18</th>
      <td>210870</td>
      <td>0.01</td>
      <td>0.01</td>
      <td>0.036147</td>
    </tr>
    <tr>
      <th>19</th>
      <td>213966</td>
      <td>0.01</td>
      <td>0.01</td>
      <td>0.036147</td>
    </tr>
    <tr>
      <th>20</th>
      <td>225119</td>
      <td>0.01</td>
      <td>0.01</td>
      <td>0.036147</td>
    </tr>
    <tr>
      <th>21</th>
      <td>241596</td>
      <td>0.01</td>
      <td>0.01</td>
      <td>0.036147</td>
    </tr>
    <tr>
      <th>22</th>
      <td>243294</td>
      <td>0.01</td>
      <td>0.01</td>
      <td>0.036147</td>
    </tr>
    <tr>
      <th>23</th>
      <td>246897</td>
      <td>0.01</td>
      <td>0.01</td>
      <td>0.036147</td>
    </tr>
    <tr>
      <th>24</th>
      <td>251658</td>
      <td>0.01</td>
      <td>0.01</td>
      <td>0.036147</td>
    </tr>
    <tr>
      <th>25</th>
      <td>25202</td>
      <td>0.01</td>
      <td>0.01</td>
      <td>0.036147</td>
    </tr>
    <tr>
      <th>26</th>
      <td>252915</td>
      <td>0.01</td>
      <td>0.01</td>
      <td>0.036147</td>
    </tr>
    <tr>
      <th>27</th>
      <td>280935</td>
      <td>0.01</td>
      <td>0.01</td>
      <td>0.036147</td>
    </tr>
    <tr>
      <th>28</th>
      <td>3</td>
      <td>0.03</td>
      <td>0.03</td>
      <td>0.072627</td>
    </tr>
    <tr>
      <th>29</th>
      <td>53625</td>
      <td>0.01</td>
      <td>0.01</td>
      <td>0.036147</td>
    </tr>
    <tr>
      <th>30</th>
      <td>54582</td>
      <td>0.01</td>
      <td>0.01</td>
      <td>0.036147</td>
    </tr>
    <tr>
      <th>31</th>
      <td>64930</td>
      <td>0.01</td>
      <td>0.01</td>
      <td>0.036147</td>
    </tr>
    <tr>
      <th>32</th>
      <td>6548</td>
      <td>0.05</td>
      <td>0.05</td>
      <td>0.127317</td>
    </tr>
    <tr>
      <th>33</th>
      <td>73764</td>
      <td>0.01</td>
      <td>0.01</td>
      <td>0.036147</td>
    </tr>
    <tr>
      <th>34</th>
      <td>74361</td>
      <td>0.03</td>
      <td>0.03</td>
      <td>0.072726</td>
    </tr>
    <tr>
      <th>35</th>
      <td>84477</td>
      <td>0.01</td>
      <td>0.01</td>
      <td>0.036147</td>
    </tr>
    <tr>
      <th>36</th>
      <td>98628</td>
      <td>0.01</td>
      <td>0.01</td>
      <td>0.036147</td>
    </tr>
  </tbody>
</table>
</div>



## Nota: 

Para obtener los mismos resultados con tres algoritmos (NetworkX, método de la potencia adaptado, y método de la potencia visto en clase), fue necesario generar la matriz de Google a partir de la matriz de adyacencias. Los algoritmos del método de la potencia se aplicaron sobre la matriz de Google para obtener su eigenvector, de esta manera, los 3 algoritmos identificaron a la página "2" como la de mayor rank.


La matriz de Google considera el factor de amortiguación a los nodos que no tienen outlinks, lo que ayuda a que no se vean beneficiados por lo mismo. Sin aplicar la matriz de Google, los métodos de la potencia identificaban al nodo "74361" como el de mayor rank; se puede apreciar en el grafo que ese nodo no cuenta con outlinks, y sí cuenta con varios inlinks. Sin embargo, intuitivamente a partir del grafo, se puede apreciar que el nodo "2" debería ser el que tiene mayor Rank. La función de `pagerank` de NetworkX ya considera un factor de amortiguación por default, por lo que arrojaba un resultado distinto al de los métodos de la potencia sin ser aplicados a la matriz de Google.


Desgraciadamente, no fue posible utilizar el algoritmo desarrollado por el equipo durante la práctia 1.2 para obtener los eigenvalores a partir de la factorización QR para la matriz de adyacencias de 37 páginas.

## Obtenemos con el algoritmo QR implementado en la Práctica 1.2


```python
eigval37, eigvec37 = Eig_Householder(A_37)
```

## Referencias

- [Google Matrix](https://networkx.org/documentation/stable//reference/algorithms/generated/networkx.algorithms.link_analysis.pagerank_alg.google_matrix.html#networkx.algorithms.link_analysis.pagerank_alg.google_matrix)
- [Código QR](https://github.com/optimizacion-2020-2-gh-classroom/practica-1-segunda-parte-graph-caroacostatovany)
- [PageRank](https://networkx.org/documentation/stable//reference/algorithms/generated/networkx.algorithms.link_analysis.pagerank_alg.pagerank.html)
- [Etiquetas para los nodos](https://stackoverflow.com/questions/14665767/networkx-specific-nodes-labeling)
- [set de datos de Stanford](http://snap.stanford.edu/data/web-Stanford.html)

06edd5bbf8d03067280c42d63de9ea7f59e7d393
ipynb
Jupyter Notebook
Implementacion/PageRank_implementacion.ipynb
lecepe00/PageRank-eigenvectores-eigenvalores
de8381b49784b1ddea058d7b4828cc709aae7c1f
MIT
Implementacion/PageRank_implementacion.ipynb
lecepe00/PageRank-eigenvectores-eigenvalores
de8381b49784b1ddea058d7b4828cc709aae7c1f
MIT
2020-11-30T04:16:10.000Z
2020-12-09T23:33:28.000Z
Implementacion/PageRank_implementacion.ipynb
lecepe00/PageRank-eigenvectores-eigenvalores
de8381b49784b1ddea058d7b4828cc709aae7c1f
MIT
2020-12-09T01:24:19.000Z
2020-12-09T01:24:19.000Z
Qwen/Qwen-72B
1. YES
2. YES