Definició i notacions

Calcular la inversa d’una matriu quadrada en R és ben fàcil, fem servir la funció solve(). Més endavant entendreu perquè calcular la inversa és el mateix que “resoldre”.

A <- matrix(c(2, 0, 3, 4, 7, 5,  5, 4, 9), ncol = 3)
A
     [,1] [,2] [,3]
[1,]    2    4    5
[2,]    0    7    4
[3,]    3    5    9
solve(A)
           [,1]        [,2]       [,3]
[1,]  1.4827586 -0.37931034 -0.6551724
[2,]  0.4137931  0.10344828 -0.2758621
[3,] -0.7241379  0.06896552  0.4827586

R fa sempre els càlculs amb la versió decimal (se’n diu de coma flotant, vegeu el final de l’anterior guia) dels nombres implicats. Entre altres coses, això té l’aventatge de ser molt més eficient. I no importa massa la dimensió de la matriu.

A <- matrix(sample(-24:24), nrow = 7) # són els nombres des de -24 fins a 24 reordenats aleatoriament
iA <- solve(A)
round(iA,4)
        [,1]    [,2]    [,3]    [,4]    [,5]    [,6]    [,7]
[1,] -0.0174  0.0123  0.0152  0.0113 -0.0051 -0.0076 -0.0199
[2,]  0.0117  0.0127  0.0222  0.0008  0.0091 -0.0138  0.0137
[3,]  0.0065 -0.0017  0.0005  0.0174 -0.0113 -0.0285 -0.0238
[4,] -0.0136 -0.0077 -0.0072  0.0136  0.0184 -0.0084  0.0133
[5,]  0.0219  0.0409 -0.0044 -0.0099 -0.0284  0.0029  0.0585
[6,] -0.0263 -0.0209  0.0034 -0.0030 -0.0095 -0.0247 -0.0257
[7,]  0.0011 -0.0214  0.0101  0.0132  0.0113 -0.0030  0.0244
A %*% iA
              [,1]          [,2]          [,3]          [,4]          [,5]
[1,]  1.000000e+00 -4.857226e-17  9.020562e-17  2.740863e-16 -1.804112e-16
[2,] -8.326673e-17  1.000000e+00  0.000000e+00  0.000000e+00  1.665335e-16
[3,]  5.204170e-17  5.551115e-17  1.000000e+00 -5.551115e-17  2.775558e-17
[4,]  7.285839e-17  0.000000e+00  0.000000e+00  1.000000e+00  8.326673e-17
[5,] -2.775558e-17  5.551115e-17 -5.551115e-17  0.000000e+00  1.000000e+00
[6,]  7.979728e-17 -5.551115e-17  0.000000e+00 -1.110223e-16  1.665335e-16
[7,]  1.040834e-17 -2.775558e-17  0.000000e+00  0.000000e+00  6.938894e-17
              [,6]          [,7]
[1,] -1.908196e-17 -8.326673e-17
[2,] -1.387779e-16 -1.110223e-16
[3,] -1.387779e-17  0.000000e+00
[4,]  2.775558e-17  0.000000e+00
[5,]  3.469447e-17 -5.551115e-17
[6,]  1.000000e+00  1.665335e-16
[7,]  4.857226e-17  1.000000e+00

Com veieu, aixó té l’inconvenient que de vegades no s’entén a la primera el què hem obtingut, aquesta matriu identitat no ho sembla gaire degut als petitissims errors que s’acumulen en els càlculs amb decimals.

Podem demanar el resultat arrodonit i ho veiem clar:

round(A %*% iA,4)
     [,1] [,2] [,3] [,4] [,5] [,6] [,7]
[1,]    1    0    0    0    0    0    0
[2,]    0    1    0    0    0    0    0
[3,]    0    0    1    0    0    0    0
[4,]    0    0    0    1    0    0    0
[5,]    0    0    0    0    1    0    0
[6,]    0    0    0    0    0    1    0
[7,]    0    0    0    0    0    0    1

Resoldre un sistema d’equacions lineals és multiplicar per la inversa.

Com que un sistema d’equacions lineals es pot escriure matricialment com a \({\mathbf A}{\mathbf X}={\mathbf b}\), la solució serà \({\mathbf X}={\mathbf A}^{-1}{\mathbf b}\). Per tant, per resoldre el sistema \[\begin{array}{lcc} x+2y-z&=&-5\\2x-y+z&=&6\\x-y-3z&=&-3 \end{array}\] farem:

A <- matrix(c(1,2,1,2,-1,-1,-1,1,-3), nrow = 3)
A
     [,1] [,2] [,3]
[1,]    1    2   -1
[2,]    2   -1    1
[3,]    1   -1   -3
b <- c(-5, 6, -3)
solve(A,b)
[1]  1 -2  2

Observa que la funció solve() que hem fet servir per calcular la inversa, si li donem un segon argument, ens torna el producte de la inversa per aquest segon argument.

solve(A) %*% b- solve(A,b)  # pràcticament zero, és el mateix
              [,1]
[1,] -2.220446e-16
[2,]  0.000000e+00
[3,]  0.000000e+00

Matrius ortogonals

La matriu de l’exemple de la guia:

P<- matrix(c(1,-1,1,1),nrow = 2) / sqrt(2)
P
           [,1]      [,2]
[1,]  0.7071068 0.7071068
[2,] -0.7071068 0.7071068
solve(P) # la inversa és la transposada
          [,1]       [,2]
[1,] 0.7071068 -0.7071068
[2,] 0.7071068  0.7071068
P %*% t(P)
     [,1] [,2]
[1,]    1    0
[2,]    0    1
t(P) %*% P
     [,1] [,2]
[1,]    1    0
[2,]    0    1

Per ortogonalitzar una matriu, R fa servir un mètode anomenat de Gram-Schmidt, disponible al paquet pracma. No hi entrem en aquest curs, però el podeu trobar força ben explicat a (https://ca.wikipedia.org/wiki/Proc%C3%A9s_d%27ortogonalitzaci%C3%B3_de_Gram-Schmidt).

L’apliquem a un exemple.

library(pracma) # potser us cal fer install.packages("pracma") si no el teniu instal·lat

A <- matrix(sample(-12:12), nrow = 5) # són els nombres des de -12 fins a 12 reordenats aleatoriament
A
     [,1] [,2] [,3] [,4] [,5]
[1,]   -5   -8    9    2    6
[2,]   10    3    0   -1   12
[3,]    5   11   -4    1  -12
[4,]   -3    4   -7   -9   -2
[5,]    8  -11   -6  -10    7
nA <- gramSchmidt(A)$Q
nA
           [,1]       [,2]        [,3]         [,4]        [,5]
[1,] -0.3348248 -0.4106521  0.47371293 -0.538244532 -0.45293062
[2,]  0.6696495  0.1037150  0.39569159 -0.478767274  0.39372930
[3,]  0.3348248  0.5762496  0.01891833 -0.004124846 -0.74528780
[4,] -0.2008949  0.2393614 -0.64384371 -0.693568338  0.08231416
[5,]  0.5357196 -0.6566968 -0.45180926  0.004546161 -0.27857058
round(nA %*% t(nA), 4)
     [,1] [,2] [,3] [,4] [,5]
[1,]    1    0    0    0    0
[2,]    0    1    0    0    0
[3,]    0    0    1    0    0
[4,]    0    0    0    1    0
[5,]    0    0    0    0    1

Els models input-output de Leontieff

Podeu consultar el llibre Matemáticas aplicadas a la Administración y a la Economía de J.C. Arya i L.W. Lardner, Pearson, México 2009. Algunes pàgines al fitxer [http://pascal.upf.edu/math21/mates1-21/matrius/Arya-Lardner-p.pdf] Arya-Lardner-p.pdf.

Les dades reals referents a Catalunya les podeu trobar a [https://www.idescat.cat/estad/mioc] però com a exercici del curs són massa complexes d’explicar i analitzar ja que considera desenes de sectors per analitzar l’economia del país. O les dades per a l’economia mundial [http://www.wiod.org/home].

Aquí ens referirem a la pàgina [https://en.wikipedia.org/wiki/Input%E2%80%93output_model] Model input-output de la Viquipedia.

Si dividim la economia d’un país en \(N\) sectors, la matriu input-output té elements \(a_{ij}\) que s’interpreten com que el sector \(j\) per produir una unitat necessita consumir \(a_{ij}\) unitats de cada sector \(i\).

Per tant, si \(\mathbf{x} = (x_1, x_2, \dots, x_N)'\) és el vector de output dels sectors, i \(\mathbf{d}=(d_1, d_2, \dots, d_N)'\) és el vector de les demandes finals de cada sector, es tindrà \[x_i = a_{i1} x_1 + a_{i2} x_2, \dots, a_{iN} x_N + d_i\] o en forma matricial \[\mathbf{x} = \mathbf{A} \mathbf{x} + \mathbf{d}\] o bé, posant \(\mathbf{x}=\mathbf{I}\mathbf{x}\) amb \(\mathbf{I}\) la matriu identitat, \[(\mathbf{I}-\mathbf{A})\mathbf{x} = \mathbf{d}, \qquad \mathbf{x}=(\mathbf{I}-\mathbf{A})^{-1} \mathbf{d}\] suposant que la matriu és invertible

Amb les dades de l’exemple de la Viquipèdia, !(input-output-viquipedia.png)

tindrem \(\mathbf{A}=\left(\begin{array}{ccc} 5/90&15/80&2/30\\10/90&20/80&10/30\\10/90&15/80&5/30\end{array}\right)\), i \(\mathbf{d}=(68,40,0)'\).

A = matrix(c(5/90, 15/80, 2/30, 10/90, 20/80, 10/30, 10/90, 15/80, 5/30), ncol=3, byrow=TRUE)
A
           [,1]   [,2]       [,3]
[1,] 0.05555556 0.1875 0.06666667
[2,] 0.11111111 0.2500 0.33333333
[3,] 0.11111111 0.1875 0.16666667

Amb la qual cosa la matriu inversa de Leontieff és

L = solve(diag(1, nrow=3)-A)
L
          [,1]      [,2]      [,3]
[1,] 1.1250000 0.3375000 0.2250000
[2,] 0.2592593 1.5592593 0.6444444
[3,] 0.2083333 0.3958333 1.3750000

que permet calcular la producció necessària per cobrir una demanda determinada, per exemple si prenem la demanda que hem introduït al model:

L%*% c(68,40,0)
     [,1]
[1,]   90
[2,]   80
[3,]   30

Però si preveiem que la demanda serà \((120,70,10)'\), per satisfer-la la producció haurà de ser:

L%*% c(120,70,10)
          [,1]
[1,] 160.87500
[2,] 146.70370
[3,]  66.45833

Exercici

Donades les dades de consum i producció d’una economia amb tres sectors S1, S2, S3, calcula la matriu input-output de manera semblant a com ho hem fet suara, calcula la matriu inversa de Leontieff i fes-la servir per obtenir els nivells de producció necessaris per satisfer una demanda de \((80, 200, 500)'\)

Tingues en compte que no s’han donat els nivells totals de producció (que han de coincidir amb els inputs totals) i també hauràs de calcular els nivells de inputs auxiliars.

Respostes: la matriu input-output és (en columnes) 0.4, 0.2, 0.2, 0.5, 0.1, 0.2, 0.3, 0.1, 0.1. La inversa de Leontieff, també en columnes 2.56, 0.65, 0.71, 1.65, 1.56, 0.71, 1.03, 0.39, 1.43 i la producció demanada és 1048.8, 555.8, 912.13.

Exercici més difícil

Donades les dades de la figura, que també pots trobar a [Exemple-IO-cat.xlsx], se’t demana construir la matriu input-output corresponent i els vectors de producció total i de demanda total.

Ajuda: La tercera fila de la matriu ha de ser 0.0071, 0.0062, 0.2497, 0.0147, 0.0039, 0.0108

Calcula també la matriu inversa de Leontieff. (Ajuda, la tercera columna ha de ser 0, 0, 1, 0, 0, 0).

Finalment, calcula els nivells de producció necessaris per satisfer una demanda de c(5000, 150000, 30000, 80000,18000, 10000). (Ajuda: el resultat ha de ser 14724.48, 257186.13, 44937.54, 119037.51, 23963.96, 16767.53).