Módulo I: (c) Soluciones de SELS (Introducción a la formulación y resolución de Modelos Matemáticos Matriciales) - leangior/MAP GitHub Wiki

Matemática Aplicada (Introducción a la Modelación Matemática en Cs. Ambientales)

  • Comisión Ecología: Dr. Leandro Giordano / Prof. Silvana Ávila

Módulo I: Sistemas de Ecuaciones Lineales y Modelos Matriciales

Objetivo

  • Aprender a formular y operar modelos matemáticos para la estimación de magnitudes físicas, ecológicas o demográficas, mediante el uso de sistemas de ecuaciones lineales y su representación matricial

Aplicaciones

  1. Hallar el conjunto soluciones posibles para un sistema de ecuaciones lineales que caracteriza una situación específica en tiempo o espacio [modelos estáticos]

    • Ejemplos: Estimar la distribución de individuos de distintas poblaciones en distintos sitios de acuerdo al stock disponible en cada sitio y a la cantidades necesarias por cada individuo de la población, interpolación espacial o temporal de valores de magnitudes físicas o ambientales formulando los correspondientes sistemas de ecuaciones lineales para determinar coeficientes de polinomios o planos (e.g. considerando la variación de estas magnitudes sobre una transecta o sobre una superficie)
  1. Formular modelos dinámicos en tiempo discreto y evaluar su operación [modelos dinámicos]

    • Ejemplos: Identificar la estructura de un sistema demográfico y formular las ecuaciones lineales que representan la dinámica del tránsito de una cohorte, identificando la matriz de transición del sistema y, de ahí, realizar proyecciones fundadas sobre la evolución del sistema (e.g. concepto de autovalor)
  1. Formular y operar sistemas discretos lineales para la representación de fenómenos de propagación o tránsito de una señal [funciones de transferencia de sistemas lineales]

    • Ejemplos: Simular la propagación del vuelco de una sustancia contaminante entre 2 puntos situados en una red de drenaje o simular el tránsito de hidrogramas

Soluciones de Sistemas de Ecuaciones Lineales

Matriz Identidad

Es sabido que el número 1 desempeña un papel fundamental en los conjuntos de números. En efecto, se sabe que el producto entre cualquier número real k y el número 1 es, siempre, el primero (k.1=k). Esto se conoce como propiedad de identidad del número 1. Luego, cabe preguntarse ¿Existirá algún tipo de matriz dentro del Álgebra de matrices que bien pueda cumplir con esta propiedad? Esto es ¿existirá algún tipo de matriz que sea algo así como el número 1 en las matrices?

En efecto, esto es así. Más específicamente, se introduce la definición de un tipo especial de matriz cuadrada, la matriz identidad de orden n.

Sea In una matriz de dimensiones nxn, cuyas entradas o elementos de la diagonal principal son iguales a 1 y el resto de los elementos es igual a 0. En notación expandida, generalizando esto es:

$$\mathbf{I_n}=\begin{align} \begin{bmatrix} 1 & 0 & . &. & . & 0 \\ 0 & 1 & . &. & . & 0 \\ . & . & . &. & . & . \\ . & . & . &. & . & . \\ . & . & . &. & . & . \\ 0 & 0 & . &. & . & 1 \end{bmatrix} \end{align}$$

diremos que In es la matriz identidad de orden n.

Así, por ejemplo podemos establecer que:

$$\begin{align} \mathbf{I_2}= \begin{bmatrix} 1 & 0 \\ 0 & 1 \end{bmatrix} \; \; \; \mathbf{I_3}= \begin{bmatrix} 1 & 0 & 0 \\ 0 & 1 & 0 \\ 0 & 0 & 1 \\ \end{bmatrix} \; \; \; \mathbf{I_5}= \begin{bmatrix} 1 & 0 & 0 & 0 & 0 \\ 0 & 1 & 0 & 0 & 0 \\ 0 & 0 & 1 & 0 & 0 \\ 0 & 0 & 0 & 1 & 0 \\ 0 & 0 & 0 & 0 & 1 \end{bmatrix} \end{align}$$

Luego, además puede observarse la propiedad señalada previamente (¿Existe lo que es un 1 a los números en el Álgebra de las matrices?). Sea A una matriz de dimensiones mxn, luego es relativamente sencillo verificar que:

$$\mathbf{I_m A}=\mathbf{A I_n}=\mathbf{A}$$

Más adelante, cuando introduzcamos el desarrollo del método de eliminación valiéndonos de las propiedades de equivalencia de matrices y SELS y aplicando la reducción por filas a la forma escalonada, veremos que esta igualdad nos brindará un soporte teórico robusto, más específicamente para la resolución de sistemas compatibles y determinados.

Ahora bien, suponga que A es una matriz cuadrada (de dimensiones nxn). Se definen las potencias de una matriz:

$$\mathbf{A^p}=A.A.A.....A \; \; (p \; veces)$$

y se establece que

$$\mathbf{A^0}=\mathbf{I_n}$$

Más adelante, veremos que en Sistemas Dinámicos Lineales la potenciación adquiere especial interés

Por ejemplo, para el Sistema Lineal $\mathbf{A.x_{(p)}=x_{(p+1)}}$, en donde A es una matriz de orden nxn (denominada de transición o transformación) y xp es un n-vector de n variables (denominadas variables de estado) que varían entre dos intervalos de tiempo p y p+1 (de acuerdo a la información contenida en A), puede demostrarse que determinados los n-valores iniciales de x(0) para cualquier valor entero p que $$\mathbf{A^p x_{(0)}=x_{(p)}}$$.

Matriz Aumentada (Ampliada)

Sea

$$\mathbf{Ax=b}$$

la representación matricial de un SEL de mxn. Se denomina matriz aumentada o ampliada [A|b], a la matriz de dimensiones mx(n+1) que se obtiene al combinar la matriz de coeficientes A de mxn con b, el m-vector columna de términos independientes, de la siguiente manera:

$$\begin{align} \mathbf{[A|b]}=\begin{bmatrix} a_{11} & a_{12} & . &. &. & a_{1n} &|& b_{11} \\ a_{21} & a_{22} & . &. &. & a_{2n} &|& b_{21} \\ . & . & . &. & . &. &|& . \\ . & . & . &. & . &. &|& . \\ . & . & . &. & . &. &|& . \\ a_{m1} & a_{m2} & . &. &. & a_{mn} &|& b_{m1} \end{bmatrix} \end{align}$$

Matriz con Forma Escalonada Reducida por Filas

A fin de sistematizar el método de eliminación de incógnitas, se propone iniciar el cómputo elaborando la matriz aumentada [A|b] (agrupando la matriz de coeficientes y el vector de términos independientes), y mediante un conjunto de operaciones sobre la base de los principios de proporcionalidad y aditividad, para luego obtener una matriz equivalente por filas [C|d], de manera tal que el SELs representado sea más sencillo de resolver.

Por ejemplo, si

$$\begin{align} \mathbf{[C|d]}=\begin{bmatrix} 1 & 1 & 0 & 2 &|& 4 \\ 0 & 1 & 0 & -1 &|& -5 \\ 0 & 0 & 1 & 3 &|& 6 \\ \end{bmatrix} \end{align}$$

representa la matriz aumentada de un sistema lineal, ya reducida a una forma simple mediante operaciones elementales basadas en los principios de linealidad, podremos observar que es relativamente sencilla la solución a partir de las ecuaciones correspondientes (podemos definir $x_4$ en función de $x_3$ y reemplazar en el resto de las ecuaciones)

$$\begin{align} \begin{matrix} x_1&+&x_2&+&0&+&2x_4&=&4 \\ 0&+&x_2&+&0&-&x_4&=&-5 \\ 0&+&0&+&x_3&+&3x_4&=&6 \end{matrix} \end{align}$$

Sol: Es un sistema compatible e indeterminado con $s=(\dfrac{2x_3}{3},-1(3+\dfrac{x_3}{3}),x_3,2-\dfrac{x_3}{3})$ [verifique]

Luego, se podría intuir que la manipulación de la matriz aumentada mediante un conjunto de reglas simples (brindadas por los principios de linealidad) podría permitirnos obtener formas simples como la precedente, de manera tal que puedan deducirse fácilmente las soluciones y más aun si el sistema es compatible y determinado.

En efecto, la forma simple que observamos en la matriz precedente se denomina forma escalonada reducida por filas, para lo cual debe satisfacer las siguientes propiedades:

  • Todas las filas que constan sólo de ceros, si las hay, estarán en la parte inferior de la matriz.
  • La primera entrada o elemento distinto de cero en una fila, al leer de izquierda a derecha, es un uno. Esta entrada o elemento se denomina uno principial de su fila
  • Para cada fila que no consta de ceros, el uno principal aparece a la derecha y abajo de cualquier otro uno principal de las filas precedentes.
  • Si una columna contiene un uno principal, el resto de las entradas o elementos de esta son iguales a cero.

Finalmente si se parte de una matriz aumentada [A|b] y se desea obtener una matriz equivalente con forma escalonada reducida por filas [C|d], sobre la base de los principios de proporcionalidad y aditividad pueden definirse las siguientes operaciones elementales por fila:

  • Intercambiar filas (reemplazar la i-ésima fila por la r-ésima fila).
  • Multiplicar una fila por un valor escalar $k$ distinto de 0 (multiplicar la i-ésima fila por un valor $k$).
  • Sumar $k$ veces la r-ésima fila a la i-ésima fila, siendo i y r distintos.

Cualquier operación elemental por fila terminará afectando los valores de los coeficientes de dicha fila. Luego si A1 es una matriz de mxn, si se intercambia la fila 1 por la 2, en la nueva matriz equivalente obtenida A2, los elementos de la primera fila serán los elementos de la segunda fila de A1 y los de la segunda los de la primera de A1. También, si se multiplica por 2 la tercera fila de A1, en la matriz equivalente resultante A2 la tercera fila estará compuesta por los elementos resultantes de este producto. Finalmente, si a la m-ésima fila de A1 se le suma k veces la i-ésima fila, la m-ésima fila de A2 estará conformada por esta suma de vectores.

La forma escalonada reducida por filas en SELs compatibles y determinados

Sea Ax=b un SEL compatible y determinado de dimensiones nxn. De acuerdo a lo previamente definido, se procede a resolver el sistema representando la matriz ampliada [A|b]. Veremos que si se procede con cautela, en este caso la reducción por filas a la forma escalonada, si se cumplen las propiedades previamente enunciadas, conduce a una matriz ampliada [C|d], en donde C será igual a In, la matriz identinad de orden n.

Luego, utilizando las reglas elementales del Álgebra de matrices que se han desarrollado, sea Ax=b un SEL compatible y determinado de dimensiones de nxn y [C|d] la forma escalonada reducida por filas de [A|b], podemos formular:

$$\begin{align} \mathbf{Cx=d} \\ \mathbf{C=I_n} \\ \end{align}$$

Lo que finalmente conduce a:

$$\begin{align} \mathbf{I_n x=d} \\ \mathbf{x=d} \\ \end{align}$$

Igualdad que nos permite concluir que el vector solución del sistema Ax=b está constituido por los elementos del vector d de términos independientes de la matriz [C|d], de forma escalonada reducida por filas (i.e. los elementos que conforman la columna n+1).

Matriz inversa y solución de SELs compatibles y determinados

Sea A una matriz cuadrada de dimensiones nxn. Se define a la inversa de la matriz A-1, como aquella que satisface la siguiente relación

$$\begin{align} \mathbf{A^{-1}A=AA^{-1}=I_n} \end{align}$$

Si una matriz presenta inversa se dice que es invertible o no singular (esto ocurre sólo para matrices cuadradas de dimensiones nxn, en donde todas las filas sean linealmente independientes).

De acuerdo a la definición precedente, se establece que la solución de un SEL determinado y compatible Ax=b (esto es, con n ecuaciones linealmente independientes y n incógnitas) puede obtenerse aplicando el producto de la matriz de la inversa de la matriz de coeficientes A en ambos términos de la igualdad, de manera tal que operando se obtiene:

$$\begin{align} \mathbf{Ax}&=\mathbf{b} \\ \mathbf{A^{-1}Ax}&=\mathbf{A^{-1}b} \\ \mathbf{I_{n}x}&=\mathbf{A^{-1}b} \end{align}$$

Lo que conduce a la solución del sistema, expresada mediante la operación del producto de A-1 y b:

$$\begin{align} \mathbf{x}&=\mathbf{A^{-1}b} \end{align}$$

En estos dos últimos puntos se ha podido ver cómo se vinculan la forma escalonada reducida por filas, la matriz identidad, la inversa de la matriz de coeficientes y el vector solución de un SEL compatible y determinado Ax=b. Esto es relevante, puesto que gracias a lo primero veremos que podremos valernos de la reducción mediante operaciones elementales por filas para obtener una matriz C=In en la forma [C|d], de tal manera que x=d es el vector solución de Ax=b. En relación a lo segundo, al menos la mayoría de las aplicaciones de planillas de cálculo (escritorio o web) permiten realizar el cómputo de la matriz inversa (ver aquí), tanto como el producto entre matrices (ver aquí), de manera tal que conociendo A y b puede hallarse fácilmente x operando en una planilla de cálculo (lo que puede resultar muy útil para sistemas compatibles y determinados de dimensiones grandes)

Formulación y Resolución de Modelos Lineales (Matriciales)

En esta sección primeramente veremos como aplicar los principios teóricos expuestos en las secciones precedentes. Ya sea para la formulación de modelos lineales mediante sistemas de ecuaciones, tanto como para su posterior resolución mediante la representación matricial (producto matricial y matriz aumentada) y la aplicación de la reducción por filas de la matriz aumentada a una forma equivalente y escalonada, a fin de hallar las posibles soluciones, generalmente en el contexto de estimación de magnitudes que caracterizan una situación específica en tiempo o espacio (modelos estáticos).

Luego, se presentarán algunos elementos teóricos para la formulación de Modelos Dinámicos Lineales (Matriciales), particularmente presentando formulaciones clásicas en el contexto de Ecología de Poblaciones mediante el uso de transformaciones lineales. Asimismo, se verá la utilidad del cómputo del determinante para la estimación de autovalores, los cuales serán útiles para la interpretación de la dinámica de estos sistemas (i.e. establecer si son estacionarios, inflacionarios, deflacionarios). Específicamente, para este tema en particular se dispone de una ficha elaborada por la Comisión, a la que referiremos en las ocasiones que sean necesarias a fin de simplificar la redundancia (ver artículo). Finalmente, se introducirá el uso de funciones de transferencia u operadores lineales, aproximando elementos de la teoría de sistemas lineales y su aplicación en procesos de propagación o tránsito de energía/masa.

Por último, podrá notarse que el desarrollo del problema siempre se dividirá en dos componentes:

  • (a) Identificación/Visualización/Conceptualización:

    • Consiste en la identificación del objetivo perseguido, la visualización de los datos e información disponible y la elaboración del modelo matemático que permitirá abordar la resolución del problema (expresado como SEL o en su represntación matricial)
  • (b) Operación/Resolución

    • Consiste en la aplicación de los principios teóricos de resolución para la determinación de las posibles soluciones al problema

El corrrecto desarrollo de ambas componentes conducirá a un resultado satisfactorio, el error operacional puede ser tanto más frecuente que el conceptual. Un error conceptual no conducirá a una resolución satisfactoria, pues la operación se apoya sobre la conceptualización (al menos en el ámbito de la modelación).

Modelos Estáticos

En el transcurso de este taller, se entenderá como modelo estático (lineal) a todo aquel que caracteriza una situación específica en tiempo o espacio. Así, un modelo estático lineal estará conformado por un SEL de m ecuaciones y n incógnitas (con mxn coeficientes y m términos independientes), utilizadas para describir las relaciones entre magnitudes de un fenómeno observado, en un intervalo tiempo o instante específico o para una situación de régimen permanente, ya sea de forma distribuida (en varios puntos) o agregada (concentrada en un punto) en el espacio de interés. En pocas palabras, sin considerar al tiempo como variable (ya sea por el caráter teórico de la situación - e.g. régimen permanente - o por que el interés esté centrado en el cómputo de magnitudes de atributos en una situación temporal específica).

Ejemplo de Aplicación 1. Modelo de Distribución de Temperatura para una placa aislada y con bordes a temperatura constante

(a) Planteo del Problema (identificación/visualización/conceptualización)

Supóngase que se dispone de una placa perfectamente aislada térmicamente tanto por arriba y por debajo. A la vez, asuma que los bordes de la placa se encuentran a temperatura constante y que el fin perseguido consiste en estimar la temperatura en 4 puntos, igualmente espaciados, como lo muestra la Fig. 1.

Fig.1 - Condiciones observadas para formular el desarrollo de un modelo matemático lineal a fin de estimar la temperatura T1,T2,T3,T4

Para poder estimar las temperaturas, veremos que bien podría formularse un sistema lineal como modelo matemático. La definición de este modelo lineal puede desarrollarse mediante 2 pasos, que implican la asunción de criterios que impongan linealidad:

  • El primero consiste en dividir 'teóricamente' la placa en 9 nueve celdas (trazando 2 líneas horizontales y 2 verticales que pasen por los puntos en los que se desea estimar la temperatura), señalando los puntos de intersección entre estas celdas (como se muestra en la Fig. 1). Primeramente, por las condiciones observadas, a los puntos situados en los bordes se les asignará la temperatura del borde correspondiente.

  • El segundo, consiste en definir un criterio para la interpolación de los valores en aquellos puntos en donde se desconozca el valor, el cual debe ser lineal. En este caso se puede elegir aplicar la regla del promedio (media aritmética) de los 4 vecinos más próximos para estimar el valor en T1,T2,T3,T4.

En esto último, de acuerdo a la definición de media aritmética (i.e. cociente entre la suma de valores y el número de elementos considerados), puede advertirse que la interpolación de cualquiera de los valores Ti (con i=1,2,3,4) procediendo con la media aritmética como criterio es en efecto una combinación lineal de los valores de los 4 vecinos más próximos y, de ahí, un método de interpolación lineal. Además, al utilizar los 4 vecinos más próximos y tener que desarrollar 4 ecuaciones para la interpolación de cada uno de los valores Ti (con i=1,2,3,4), queda claro que esto da lugar a un SEL de 4x4.

Por ejemplo, en principio se postulan las ecuaciones que definen la estimación en cada punto, en términos del valor medio de los 4 puntos más próximos. Así, de acuerdo a la información disponible (métodos y datos) podemos establecer que las estimaciones serán:

$$\begin{align} T_1&=\dfrac{60+100+T_2+T_3}{4} \\ T_2&=\dfrac{T_1+100+40+T_4}{4} \\ T_3&=\dfrac{60+T_1+T_4+0}{4} \\ T_4&=\dfrac{T_3+T_2+40+0}{4} \end{align}$$

Esto es, un sistema de ecuaciones lineales de 4x4. Reordenando, esto resulta aun más evidente:

$$\begin{matrix} 4T_1&-&T_2&-&T_3&+&0&=160 \\ -T_1&+&4T_2&+&0&-&T_4&=140 \\ -T_1&+&0&+&4T_3&-&T_4&=60 \\ 0&-&T_2&-&T_3&+&4T_4&=40 \end{matrix}$$

En consecuencia, este SEL de 4x4 es el modelo matemático que utilizaremos para estimar Ti (con i=1,2,3,4), el cual se obtiene al tomar en cuenta las consideraciones previamente formuladas (i.e. condiciones de borde conocidas, adopción del criterio de la media).

(b) Resolución del problema (operación)

Luego, aplicando los principios de operación expuestos durante el desarrollo de este taller, queda claro que:

  • (a) Es un sistema compatible y determinado de 4x4 (A.T=b, siendo T el 4-vector columna de incógnitas Ti, con i=1,2,3,4)

  • (b) En consecuencia, el conjunto solución estará compuesto por una solución única y sus elementos serán las estimaciones Ti (con i=1,2,3,4), que constituyen los valores que se desean hallar.

  • (c) Luego puede procedese a representar la matriz aumentada o ampliada [A|b] y a la posterior reducción por filas hasta obtener la forma escalonada [C|d], de manera tal de reducir C hasta que sea I4 y establecer el valor de las estimaciones mediante la igualdad T=d .

Así, a continuación se representa la matriz aumentada del sistema A.T=b (verifique):

$$\begin{align} \mathbf{[A|b]}=\begin{bmatrix} 4 & -1 & -1 & 0 &|& 160 \\ -1 & 4 & 0 & -1 &|& 140 \\ -1 & 0 & 4 & -1 &|& 60 \\ 0 & -1 & -1 & 4 &|& 40 \\ \end{bmatrix} \end{align}$$

Posteriormente, se procede a realizar las operaciones elementales por fila para reducirla hasta la forma escalonada.

En primer lugar se nota una diagonal principal en la matriz de coeficientes A incluida en [A|b] (a partir de ahora se referirá a esta como diagonal principal), compuesta por todas las entradas con i=j (mismo subíndice para fila y columna), con valor 4 para todos sus elementos (de forma tal que sería fácil convertir a estos elementos en 1 si se multiplica cada fila por 1/4).

Utilizaremos la siguiente táctica para la reducción: primeramente se intentará eliminar todos los coeficientes situados por debajo de la diagonal principal y luego todos aquellos situdaos por encima de esta, valiéndonos de operaciones elementales por fila. Asimismo, siempre que sea posible se intentará reducir los elementos de la diagonal principal para que constituyan unos principales (que como se mostró anteriormente generalmente implica la operación de multiplicación de la fila por un escalar).

Así, empezamos procediendo de forma tal que los elementos por debajo del 4 de la primer fila se eliminen ('todos los elementos que se sitúen debajo de un 1 principal en una columna deben valer 0'):

  • Para esto podemos empezar multiplicando por 4 la fila 2 y sumandole la fila 1, a fin de 'actualizar' los valores de esta fila, eliminando el -1 que se sitúa debajo del 4. Asimismo, por ídem razón, multiplicamos por 4 la fila 3 y le sumamos la fila 1, para obtener una 'actualización' de la fila 3, con el primer elemento igual a 0. Así, la matriz [A|b] se reduce con estas operaciones elementales a (verifique):

$$\begin{bmatrix} 4 & -1 & -1 & 0 &|& 160 \\ 0 & 15 & -1 & -4 &|& 720 \\ 0 & -1 & 15 & -4 &|& 400 \\ 0 & -1 & -1 & 4 &|& 40 \end{bmatrix}$$

  • Luego, podemos continuar eliminando los elementos que se encuentran por debajo del 15 en la fila 2 y columna 2, a fin de simplificar la columna 2, siguiendo los principios de la forma escalonada. Para esto multiplicamos por 15 a la fila 3 y a la fila 4 y les sumamos respectivamente a cada una la fila 2, a fin de 'actualizarlas'. Así, se continua reduciendo la matriz ampliada original obteniéndose con estas operaciones elementales a (verifique):

$$\begin{bmatrix} 4 & -1 & -1 & 0 &|& 160 \\ 0 & 15 & -1 & -4 &|& 720 \\ 0 & 0 & 224 & -64 &|& 6720 \\ 0 & 0 & -16 & 56 &|& 1320 \end{bmatrix}$$

  • En este punto y prosiguiendo con el procedimiento de eliminar todos los coeficientes que se encuentren por debajo, en la columna respectiva, de cada elemento de la diagonal principal, se puede multiplicar por 14 la fila 4 y sumarle la fila 3, a fin de eliminar la entrada o elemento (4,3) que es igual a -16 (ya que -16.14=-224). Al mismo tiempo, el resultado puede divirse por 720, a fin de que la entrada (4,4) de la diagonal principal adquiera el valor 1, de forma tal de obtener con estas operaciones elementales la forma equivalente por filas (verifique):

$$\begin{bmatrix} 4 & -1 & -1 & 0 &|& 160 \\ 0 & 15 & -1 & -4 &|& 720 \\ 0 & 0 & 224 & -64 &|& 6720 \\ 0 & 0 & 0 & 1 &|& 35 \end{bmatrix}$$

Nótese que la fila 4 se encuentra completamente reducida, de acuerdo a la forma escalonada reducida por filas. Luego, podemos suponer que si los cálculos son correctos, la estimación para T4 será 35° o lo que es lo mismo T4=35°.

  • Ahora bien, llegado este punto conviene realizar el procedimiento en sentido contrario al que se venía desarrollando (ya que la última fila está completamente reducida). Esto es, se intentará eliminar los coeficientes que se encuentren por encima de la diagonal principal y establecer en 1 el valor de los elementos de esta (i.e. transformarlos en unos principales). Así, se puede sumar la fila 3 al producto entre 64 y la fila 4 (anulando el valor -64 de la fila 4) y finalmente divirlo por 224 para que la entrada (3,3) adquiera el valor 1. Por otro lado, si a la fila 2 se le suma 4 veces la fila 4 puede anularse la entrada (2,4) que se ubica encima del uno principal de dicha columna. Así se obtiene una forma aun más simplificada con estas operaciones elementales (verifique):

$$\begin{bmatrix} 4 & -1 & -1 & 0 &|& 160 \\ 0 & 15 & -1 & 0 &|& 860 \\ 0 & 0 & 1 & 0 &|& 40 \\ 0 & 0 & 0 & 1 &|& 35 \end{bmatrix}$$

Nótese que la fila 3 está completamente reducida, por lo que si los cálculos son correctos T3=40°.

  • A continuación, a fin de eliminar el resto de los coeficientes que se encuentran por encima de la diagonal principal y transformar en uno principal a los elementos restantes en esta (para obtener una submatriz C=I4 en la forma escalonada reducida por filas [C|d] de [A|b]), se propone sumarle a la fila 2 la fila 3 y dividir el resultado de esta suma por 15 (para reducir completamente la fila 2) y a la fila 1 sumarle la fila 3, a fin de anular la entrada (1,3). Así se obtiene, con estas operaciones elementales (verifique):

$$\begin{bmatrix} 4 & -1 & 0 & 0 &|& 200 \\ 0 & 1 & 0 & 0 &|& 60 \\ 0 & 0 & 1 & 0 &|& 40 \\ 0 & 0 & 0 & 1 &|& 35 \end{bmatrix}$$

Nótese que la fila 2 está completamente reducida, de forma tal que la estimación es T2=60°

  • Finalmente, siguiendo el principio de reducción adoptado queda claro que si a la fila 1 le sumamos la fila 2 y al resultado lo dividimos por 4 se anula la entrada (1,2) y la entrada (1,1) se transforma en uno principal, obteniéndose la forma escalonada reducida por filas [C|d] de [A|b] (verifique):

$$\begin{bmatrix} 1 & 0 & 0 & 0 &|& 65 \\ 0 & 1 & 0 & 0 &|& 60 \\ 0 & 0 & 1 & 0 &|& 40 \\ 0 & 0 & 0 & 1 &|& 35 \end{bmatrix}$$

La fila 1 se redujo completamente, luego la estimación es T1=65°

  • Así, aplicando la igualdad t=d ya que C=I4, el conjunto solución queda dado por:

$$s=(T_1=65°,T_2=60°,T_3=40°,T_4=35°)$$

Ejemplo de Aplicación 2. Interpolación Polinomial de una curva y(x) para estimar otra magnitud p(y)

Planteo general

Este tipo de aplicación presenta al siguiente forma:

  • En primer lugar supóngase que se conocen n pares ordenados (x1,y1),(x2,y2),...,(xn,yn) de puntos en un sistemas de coordenadas planas. En el ámbito de aplicación, estos puntos referirán a pares de valores de una variable explicativa x (e.g. distancia a un punto de referencia sobre una trayectoria, duración o permanencia) y una variable respuesta y (e.g. elevación del terreno, concentración de una sustancia, densidad de población), por ejemplo obtenidos mediante relevamiento por muestreo.

  • No se sabe exactamente cuál es la forma explícita de la función y(x) (la función matemática que vincula y con x).

  • Se desea saber el valor de y(x) para un valor de x situado entre los valores mínimos y máximos de los valores x observados en los n pares ordenados conocidos, pero para el cual no se dipone de par y conocido (problema de interpolación, si se deseara conocer el valor de f(x) para todo valor x menor a este mínimo o mayor a este máximo se estaría extrapolando y esto no es necesariamente buena práctica científica, si bien muchas veces es la única alternativa). Por conveniencia, llamemos a este valor xp.

  • Y, mas aun, se desea evaluar el valor de función p(y) conocida, que pone en relación una magnitud p con la variable respuesta y, para el valor de x=xp.

Luego, una solución propuesta para este tipo de problema consiste en elaborar un modelo matemático que permita aproximar $y(x)$. Así posteriormente se puede evaluar $p(y(x))$ para $x=x_p$, a fin de proveer una estimación. Esto último se deduce por composición de funciones, ya que se establece una relación explícita $y(x)$ y además se conoce explícitamente $p(y)$, entonces es posible obtener mediante $p(y(x))=p(x)$, p en función de x.

En efecto, la llave del problema es la determinación de una forma explícita y(x) que ajuste a los pares ordenados (x1,y1),(x2,y2),...,(xn,yn). Para esto existen distintos métodos (estadísticos, deteminísticos). Aquí se desarrollará un método determinístico que se apoya en la propiedad que establece que dados n puntos definidos en el plano con coordenadas x,y existe uno y sólo un polinomio y(x) de grado n-1 cuya trayectoria incluye a todos estos puntos. Este polinomio será utilizado para interpolar valores de y a partir de situaciones en x. Por tanto se lo denomina polinomio de interpolación.

Todo polinomio y(x) de grado n-1 está definido por la ecuación:

$$y(x)=a_{n-1}x^{n-1}+a_{n-2}x^{n-2}+...+a_1x+a_0$$

En el tipo de problema que se enfrenta se conocen los pares (xi,yi) para i=1,2,...,n. Luego, las incógnitas son los coeficientes del polinomio de interpolación an-1,an-2,...,a1,a0.

En consecuencia, si se conocen n pares ordenados (x,y) es posible desarrollar el siguiente sistemas de ecuaciones (verifique):

$$\begin{align} \begin{matrix} a_{n-1}x_1^{n-1}&+&a_{n-2}x_1^{n-2}&+&.&.&.&+&a_{1}x_1&+&a_0&=&y_1 \\ a_{n-1}x_2^{n-1}&+&a_{n-2}x_2^{n-2}&+&.&.&.&+&a_{1}x_2&+&a_0&=&y_2 \\ .&+&.&+&.&.&.&+&.&+&.&=&. \\ .&+&.&+&.&.&.&+&.&+&.&=&. \\ .&+&.&+&.&.&.&+&.&+&.&=&. \\ a_{n-1}x_n^{n-1}&+&a_{n-2}x_n^{n-2}&+&.&.&.&+&a_{1}x_n&+&a_0&=&y_n \\ \end{matrix} \end{align}$$

El cual es un sistema compatible y determinado de n ecuaciones (linealmente independientes, puesto que todos los pares x,y son distintos) y n incógnitas (los coeficientes del polinomio de interpolación). Esto último demuestra lo previamente enunciado ('existe un y sólo un polinomio...'). A la vez, el sistema puede expresarse mediante el producto matricial:

$$\begin{align} \begin{bmatrix} x_1^{n-1}& x_1^{n-2}&.&.&.&x_1&1 \\ x_2^{n-1}& x_2^{n-2}&.&.&.&x_2&1 \\ .&.&.&.&.&.&. \\ .&.&.&.&.&.&. \\ .&.&.&.&.&.&. \\ x_n^{n-1}& x_n^{n-2}&.&.&.&x_n&1 \\ \end{bmatrix} \begin{bmatrix} a_{n-1} \\ a_{n-2} \\ . \\ . \\ . \\ a_{0} \end{bmatrix}= \begin{bmatrix} y_{1} \\ y_{2} \\ . \\ . \\ . \\ y_{n} \end{bmatrix} \end{align}$$

o en notación simplificada:

$$\mathbf{Xa=y}$$

En donde a es el vector de incógnitas del SEL (coeficientes del polinomio de interpolación).

En conclusión, para la resolución del problema se procede a representar la matriz aumentada [X|y] y se la reduce mediante operaciones elementales por fila a la forma escalonada [C|d], de tal forma que es posible determinar el vector solucion mediante la igualdad a=d, obteniéndose los valores de los coeficientes del polinomio de interpolación y(x) que enlaza los pares (yi,yi) con i=1,2,...,n.

(a) Procedimiento en un caso típico (identificación/conceptualización)
  • Supóngase en una reserva natural se ha realizado un relevamiento topográfico en tres puntos a lo largo de una transecta, como muestra la Fig. 2. En abcisas, los valores x representan la distancia al origen de la transecta. En ordenadas, los valores y representan la elevación del terreno.

Fig.2 - Pares de distancia y elevación, relevados sobre la transecta. Nótese el dominio de interpolación, definido para la variable x en el intervalo [1,3] (valores mínimos y máximos de distancia en los pares disponibles, a fines de visualización las rectas perpendiculares están ligeramente corridas). Asimismo se grafica en rojo el polinomio de interpolación y(x) de grado 2 que enlaza a los 3 puntos y para el cual pretenden obtenerse los valores de los coeficientes

Luego, se sabe que la elevación para x=1 [km] es y=3 [metros], que x=2 [km] es y=4 [metros] y x=3 [km] es y=7 [metros]. O lo que es lo mismo, se conocen los pares (1,3),(2,4) y (3,7).

  • Por otro lado se sabe que la densidad de una especie arbórea, al menos en esta reserva natural, depende linealmente de la elevación, de acuerdo a la relación:

$$p_(y)=-0.02y+2$$ en donde p es la densidad (expresada en individuos/m²) e y es la elevación (en metros).

  • Finalmente, se desea estimar la densidad de esta especie par un valor x=1.5 km

Así, si n=3, se sabe que puede ajustarse un polinomio de interpolación de grado 2 (n-1). Para esto, primeramente se debe formular el SEL de 3x3 que describe tal situación (verifique):

$$\begin{align} \begin{matrix} a_{2}1&+&a_{1}1&+&a_0&=&3 \\ a_{2}4&+&a_{1}2&+&a_0&=&4 \\ a_{2}9&+&a_{1}3&+&a_0&=&7 \end{matrix} \end{align}$$

a continuación, se formula el modelo Matricial:

$$\begin{align} \begin{bmatrix} 1 & 1 & 1 \\ 4 & 2 & 1 \\ 9 & 3 & 1 \end{bmatrix} \begin{bmatrix} a_2 \\ a_1 \\ a_0 \end{bmatrix}= \begin{bmatrix} 3 \\ 4 \\ 7 \end{bmatrix} \end{align}$$

obteniéndose la matriz aumentada [X|y]:

$$\begin{align} \begin{bmatrix} 1 & 1 & 1 | 3 \\ 4 & 2 & 1 | 4\\ 9 & 3 & 1 | 7 \end{bmatrix} \end{align}$$

la cual debe reducirse por filas a la forma escalonada [C|d] a fin de obtener el vector solución para el valor de incógnitas (coeficienets del polinomio). Luego, concociendo la estructura específica de y(x) se reemplaza la misma se obtiene p(y)=p((y(x))) (reemplazando y por y(x) en la expresión p(y)) y se evalúa el resultado para x=1.5 km.

(b) Resolución (operación)

Puede operarse siguiendo la táctica del ejemplo anterior (eliminando los elementos por debajo y por encima de la diagonal principal de X por operaciones elementales por fila), de tal forma se obtiene en una serie reducida de pasos la matriz [C|d] (verifique):

$$\begin{align} \begin{bmatrix} 1 & 0 & 0 &|& 1 \\ 0 & 1 & 0 &|& -2\\ 0 & 0 & 1 &|& 4 \end{bmatrix} \end{align}$$

Así, el polinomio de interpolación para este caso es:

$$y(x)=x^2-2x+4$$

Luego, reemplazando en:

$$p(y)=-0.02y+2$$

por composición de funciones, se obtiene:

$$\begin{align} p(y(x))&=-0.02y(x)+2 \\ p(y(x))&=-0.02(x^2-2x+4)+2\\ \end{align} $$

o lo que es lo mismo:

$$\begin{align} p(x)&=-0.02(x^2-2x+4)+2 \end{align}$$

En consecuencia, se evalúa p(x) para x=1.5 km, obteniéndose:

$p(x)=-0.02y(x)+2=-0.02(1.5^2-21.5+4)+2=-0.02*3.25+2=1.935$

Así, finalmente puede afimarse que la densidad de la especie arbórea a una distancia de 1.5 km del inicio de la transecta se estima en 1.935 individuos/m².

Ejemplo de Aplicación 3. Interpolación Lineal de un campo escalar utilizando la ecuación del plano

Generalidades

Se entiende por campo escalar a toda distribución espacial de una magnitud escalar (un número), asociando un valor con coordenada z a cada punto de un espacio 2-Dimensional (2D) con coordenadas x,y. Esto es una aplicación de $R^2 \to R$, en donde para cada par x,y definido dentro de un dominio espacial de una proyección plana, se obtiene una valor z, la cota asociada del campo escalar.

Ejemplos de campo escalar son la topografía del terreno, la temperatura, la precipitación (específicamente durante el desarrollo de un evento que cubra todo el dominio espacial considerado si se quiere garantizar la hipótesis de continuidad) o magnitudes tales como algunas densidades. Por regla general, todo campo escalar es continuo en el dominio espacial para el cual se define. Luego, la técnica que se desarrollará sólo es aplicable si la variable cuyo valor se ha de interpolar puede definirse de manera continua en el dominio espacial que se ha considerado.

La Fig. 3 muestra una representación típica en 2D de un campo escalar, correspondiente a la topografía de Cuenca del Plata (de acuerdo a un Modelo Digital de Elevación del Terreno, valores en metros). En este tipo de repŕesentación a cada situación en una proyección plana horizontal con coordenadas x,y (en este caso en unidades de longitud y latitud) se le asigna un valor z de elevación del terreno (generalmente en metros), representado en una escala de colores o tonos de un color. Otra forma típica y muy útil de representar un campo escalar es mediante isolíneas o líneas que unen aquellos puntos con coordenadas x,y de igual valor para la magnitud escalar z que define el campo (posiblemente hayan utilizado alguna vez isohipsas o curvas de nivel para la determinación de cuencas hidrográficas o de la dirección del flujo superficial).

Fig.3 - Representación típica de un campo escalar mediante un mapa de colores (Modelo Digital de Elevación del Terreno para Dominio de Cuenca del Plata). Valores de la escala de colores expresados en metros sobre nivel del mar

(a) Planteo del problema (visualización/identificación/conceptualización)
  • Supóngase que se ha medido la precipitación acumulada mensual en 6 puntos de una región de interés, como se muestra en la Fig. 4. En donde los valores x,y corresponden a las coordenadas horizontales (en km) y el valor z a la precipitación observada en dicha situación geográfica (en mm).

  • A la vez, para definir el dominio de interpolación supóngase que se grafican una serie de triángulos tal cual lo muestra la Fig. 4, construyendo una red triangular irregular.

Fig.4 - Puntos de medición (nodos) y valores obtenidos para la precipitación mensual (en mm). En azul se han trazado distintos triángulos para delimitar subdominios de interpolación

  • Ahora suponga que se desea estimar la precipitación mensual para un punto cualquiera de coordenadas x,y que se sitúe dentro de alguno de estos triángulos, utilizando un método de interpolación lineal (de manera tal que el modelo adoptado para la interpolación sea un SEL)

Pues bien, en principio se puede recordar que toda ecuación de un plano tiene la forma:

$$\begin{align} a_2 x+a_1 y+a_0=z \end{align}$$

en consecuencia, la ecuación del plano es una ecuación lineal (i.e. combinación lineal, algo que ya se vio en la página a de este módulo, en esta wiki).

Así, se verá que es posible desarrollar un método de interpolación lineal basado en la ecuación de un plano. En efecto, dados 3 puntos en un espacio 3D, existirá uno y sólo un plano que los contenga (Fig. 5). En consecuencia, dados 3 puntos podrán determinarse los coeficientes a2,a1,a0 del plano que los contiene. Así se utilizará la información disponible en cada una de las aristas de cada triángulo (visualice que en todo caso cada triángulo tendrá una inclinación igual a la del plano que incluye estos 3 puntos, siendo esta la base de la interpolación lineal que se pretende desarrollar). En efecto, este método implica hallar los coeficientes de la ecuación del plano definido por cada conjunto de 3 puntos, que a la vez definen cada triángulo. En suma, cada triángulo será un subdominio de interpolación espacial con su ecuación de plano correspondiente y los coeficientes de este último estarán determinados por un SEL de 3x3.

Fig.5 - Dados 3 puntos en el espacio 3D, con coordenadas $x_i,y_i,z_i$ conocidas, es posible determinar la ecuación del único plano que los contiene a todos. Asimismo, queda claro que 3 puntos definen un triángulo sobre una superficie 2D. De ahí, que la inclinación del triángulo sea igual a la del plano que contiene a estos puntos

Inicialmente considere los 3 puntos que definen un triángulo. Asuma que (xi,yi,zi) son los datos que se poseen para cada punto en el espacio que conforman un triángulo, con subíndice i=1,2,3 (Fig. 5). Esto es, para cada arista de un triángulo con coordenadas planas xi,yi (con i=1,2,3) y con valor de precipitación zi (con i=1,2,3). A continuación podrá visualizarse que para cada uno de estos triángulos puede formularse el siguiente SEL de 3x3, siguiendo estos principios:

$$\begin{align} \begin{matrix} a_{2}x_1&+&a_{1}y_1&+&a_0&=&z_1 \\ a_{2}x_2&+&a_{1}y_2&+&a_0&=&z_2 \\ a_{2}x_3&+&a_{1}y_3&+&a_0&=&z_3 \end{matrix} \end{align}$$

Luego, su representación matricial será:

$$\begin{align} \begin{bmatrix} x_1&y_1&1 \\ x_2&y_2&1 \\ x_3&y_3&1 \\ \end{bmatrix} \begin{bmatrix} a_2 \\ a_1 \\ a_0 \end{bmatrix}= \begin{bmatrix} z_{1} \\ z_{2} \\ z_{3} \end{bmatrix} \end{align}$$

que en notación simplificada será,

$$\begin{align} \mathbf{Xa=z} \end{align}$$

en donde el vector de incógnitas a está constituído por los coeficientes de la ecuación del plano que incluye a estos 3 puntos, los términos independientes serán los valores z de los mismos y la matriz X estará compuesta por filas cuyos dos primeros elementos sean las coordenadas x,y de cada arista y el tercero asume siempre el valor 1.

Luego, como ya se mostró, se procederá a formular la matriz aumentada que denominaremos [X|z]:

$$\begin{align} \begin{bmatrix} x_1 & y_1 & 1 &|& z_1 \\ x_2 & y_2 & 1 &|& z_2\\ x_3 & y_3 & 1 &|& z_3 \end{bmatrix} \end{align}$$

para reducirla por operaciones elementales por fila a su forma escalonada [C|d], de manera tal que el vector solución será a=d, tal como se vio para el caso de interpolación polinomial.

(b) Resolución (operación)

Supóngase que se desea estimar la precipitación mensual de acuerdo a los datos de la Fig. 4 y utilizando el método propuesto para un punto situado en (x,y)=(2.5,4). Pues bien, en principio se advierte que esta localización queda dentro del triángulo definido por los nodos 1,4 y 2.

Por tanto, se procede a elaborar el SEL de 3x3 con la información provista para estos nodos y se obtiene la siguiente matriz aumentada [X|z] (verifique):

$$\begin{align} \begin{bmatrix} 3 & 5 & 1 &|& 200 \\ 2 & 4 & 1 &|& 230\\ 3 & 3 & 1 &|& 245 \end{bmatrix} \end{align}$$

Operando de acuerdo a los principios previamente señalados se obtiene la forma escalonada reducida por filas [C|d] (verifique):

$$\begin{align} \begin{bmatrix} 1 & 0 & 0 &|& -\dfrac{15}{2} \\ 0 & 1 & 0 &|& -\dfrac{45}{2}\\ 0 & 0 & 1 &|& 335 \end{bmatrix} \end{align}$$

Por lo que la ecuación del plano que incluye a estos 3 puntos será:

$$\begin{align} -\dfrac{15}{2}x-\dfrac{45}{2}y+335=z \end{align}$$

Reemplazando por los valores x=2.5 e y=4 en la ecuación precedente (con mínimas manipulaciones algebraicas), se obtiene la solución:

$$z=-\dfrac{75}{4}-90+335$$

de manera tal que

$$z=226+\dfrac{1}{4}$$

Así, finalmente se puede afirmar que el valor estimado de precipitación mensual mediante este método para el punto con coordenadas horizontales (2.5,4) será igual a 226.25 mm.

Consideraciones de implementación (aplicaciones prácticas reales)

Ciertamente este método posee gran aplicación para la interpolación puntual de valores de un campo (en un sitio y para un intervalo de tiempo o instante específico), la generación o relleno de series temporales puntuales (e.g. temperatura, precipitación, humedad del suelo) o campos sintéticos/simulados (e.g. modelos digitales de elevación, estimaciones cuantitativas de precipitación). Al respecto de su uso para la generación de modelos digitales de elevación se recomienda este breve artículo (del cual se ha adaptado la Fig. 4 al problema del ejemplo precedente!). Particularmente en el punto 4 de este artículo encontrarán aspectos fuertemente vinculados al desarrollo teórico precedente, desde el punto de vista de una disciplina fuertemente vinculada a la Ecología. Asimismo, se podrá apreciar que los principios de interpolación mediante el uso de redes triangulares irregulares se apoyan en SELs y, de ahí, que sean útiles para cualquier atributo físico que pueda considerarse continuo en el dominio de interpolación. Otro tema es cuál debiera ser la escala espacial o la densidad de nodos (cantidad de triángulos) en un dominio de interpolación como para que la hipótesis de linealidad pueda sostenerse. La intuición indica que mientras mayor sea la densidad más sostenible es la hipótesis de linealidad.

En el ejemplo precedente se calculó la precipitación mensual para un punto situado en un subdominio de interpolación específico (i.e. el triángulo conformado por los nodos 1,4 y 2 de la Fig. 4), explicando como proceder para el primer tipo de caso. Si se deseara computar el valor de precipitación para otra situación geográfica que correspondiera a otro subdominio, para ese mes, habría que repetir el procedimiento considerando los datos x,y,z observados en los nodos (aristas) que lo definen.

Asimismo, es sabido que la precipitación mensual por lo general varía a lo largo del año, por lo que es un atributo dinámico (también concentraciones de contaminantes, la humedad del suelo, el almacenamiento subterráneo de agua, entre varios). Luego, para atributos dinámicos el proceso deberá repetirse si se desea interpolar series temporales z(t) para situaciones geográficas específicas x,y, actualizando los coeficientes conforme varíen los valores z(t) en los nodos considerados (que también pueden variar! pues las redes de medición no son necesariamente permanentes...). Para atributos estáticos, como la elevación topográfica (si bien varía lo hace escalas temporales menos apreciables), el cómputo de los coeficientes de cada plano de interpolación para cada subdominio definido se realiza una sola vez, ya que para cada situación geográfica x,y (incluida en el dominio de interpolación), el valor z puede considerarse constante. Para esto existe software que permite fácilmente codificar y resolver los distintos SELs correspondientes a cada subdominio de interpolación, en todo caso debiendo definirse la situación geográfica de cada nodo (x,y) y el valor de campo escalar asociado (z), que constituyen la información de entrada al análisis.

Ejemplo de Aplicación 4. Modelos Lineales en relaciones Insumos(Recursos)-Productos(Agentes)

Supongamos que para una situación específica se requieren K insumos de un tipo específico por cada producto elaborado. Luego si se elaboran P productos es casi evidente que el cómputo del total de los insumos necesarios se realiza sobre la base del siguiente razonamiento:

K insumos / producto elaborado x P productos elaborados = KxP insumos totales (utilizados)

sin duda alguna, una relación lineal.

En relaciones de producción fabril o agropecuaria, para las cuales se conozcan las cantidades de insumos necesarios por producto, está claro que la linealidad se sostiene (ver ejemplo inicial). Por otro lado, en relaciones ecológicas este modelo lineal tiene validez en situaciones que puedan obviarse los efectos de interacción entre recursos y agentes o bien represente una situación de régimen de equilibrio, de tal manera que pueda sostenerse la hipótesis de linealidad (situaciones para las que se acepten las propiedades de proporcionalidad y aditividad).

Por ejemplo, resulta evidente que si cada unidad de pizza requiere 1/2 kg de harina y 300 cc de agua y se disponen de 1 kg y 600 cc de agua, se obtendrán 2 pizzas. Asímismo, si cada individuo de 1 especie herbívora requiere de N kg/mes de la herbácea A y M kilogramos/mes de la herbácea B para subsistir durante un mes, es evidente que K individuos necesitarán KN kg/mes de la herbácea A y KM kilogramos/mes de la herbácea B. En estas dos situaciones se sostiene la relación de escalamiento, luego puede sostenerse la linealidad.

(a) Planteo del problema y dominio de aplicación (identificación/visualización/conceptualización)

Ahora, supongamos que se conoce la relación de cantidad de insumos/producto o recursos/agente para un conjunto de sitios en el espacio, para algunas o varias situaciones en el tiempo o para distintas condiciones necesarias que se deben cumplir en un sitio específico (e.g. esto último para el caso de distintos productos en un mismo sitio y mismo intervalo de tiempo, como en el ejemplo inicial) y además se sabe cuál es la cantidad total disponible de cada insumo. Luego, procedamos a elaborar un modelo lineal sobre la base de este postulado y que vincule esta información con el total de productos/agentes elaborados o presentes en cada sitio/situación temporal o bajo las distintas condiciones necesarias establecidas. En principio, esta última información será la incógnita.

Por empezar, si se conoce la primer relación en distintos sitios/situaciones temporales/condiciones, esta información bien podría organizarse o resumirse en una matriz de insumos/producto o recursos/agente de dimensiones mxn (m tipos de insumos/productos o recursos/agente x n sitios/situaciones temporales/condiciones):

$$\begin{align} \mathbf{A}= \begin{bmatrix} a_{11} & a_{12} & . & . & . & a_{1n} \\ a_{21} & a_{22} & . & . & . & a_{2n} \\ a_{31} & a_{32} & . & . & . & a_{3n} \\ . & . & . & . & . & . \\ . & . & . & . & . & . \\ . & . & . & . & . & . \\ a_{m1} & a_{m2} & . & . & . & a_{mn} \end{bmatrix} \end{align}$$

en donde aij, con i=1,2,...,m y j=1,2,...,n es la relación insumos/producto o recursos/agente para el i-ésimo insumo/recurso y para el j-ésimo sitio (situación temporal o condición necesaria). Asimismo, la información sobre los insumos/recursos totales disponibles bien se puede organizar en un m-vector:

$$\begin{align} \mathbf{b}= \begin{bmatrix} b_{11} \\ b_{21} \\ . \\ . \\ . \\ b_{m1} \end{bmatrix} \end{align}$$

en donde bi1, con i=1,2,...,m, se corresponde con la cantidad total disponible del i-ésimo insumo/recurso. Finalmente, si llamamos x al n-vector que contiene la cantidad total de productos elaborados/agentes presentes en cada sitio/situación temporal/condición:

$$\begin{align} \mathbf{x}= \begin{bmatrix} x_{11} \\ x_{21} \\ . \\ . \\ . \\ x_{n1} \end{bmatrix} \end{align}$$

en donde xi1, con i=1,2,...,n, se corresponde con el total de productos/agentes elaborados o presentes. Dado que esta es nuestra incógnita inicial, bien puede apreciarse que combinando estos elementos se obtiene el siguiente modelo lineal:

$$\begin{align} \begin{bmatrix} a_{11} & a_{12} & . & . & . & a_{1n} \\ a_{21} & a_{22} & . & . & . & a_{2n} \\ a_{31} & a_{32} & . & . & . & a_{3n} \\ . & . & . & . & . & . \\ . & . & . & . & . & . \\ . & . & . & . & . & . \\ a_{m1} & a_{m2} & . & . & . & a_{mn} \end{bmatrix} \begin{bmatrix} x_{11} \\ x_{21} \\ . \\ . \\ . \\ x_{n1} \end{bmatrix}=\begin{bmatrix} b_{11} \\ b_{21} \\ . \\ . \\ . \\ b_{m1} \end{bmatrix} \end{align}$$

En el cual la matriz de coeficientes A es la matriz de insumos/recursos necesarios por producto/agente en cada sitio/situación temporal/condición y el vector b de términos independientes está constituido por las cantidades totales disponibles de cada insumo o recurso. Asimismo, x es el vector de incógnitas, de modo tal que la solución del sistema estará constituida por todas las combinaciones posibles de cantidades producidas o agentes presentes en cada sitio/situación temporal/condición que puedan satisfacer las relaciones lineales del sistema, esto es que determinen un consumo total de cada insumo/recurso exactamente igual al stock disponible (observado).

Consecuentemente, en caso que m=n y los vectores fila sean linealmente independientes, el SEL será compatible y determinado con una única solución: la cantidad de productos o agentes presentes en cada sitio/situación temporal/condición que dadas las relaciones insumos/recursos por producto/agente exhiban un consumo total de cada insumo/recurso exactamente igual al stock disponible (observado).

(b) Resolución (operación)

Ciertamente, como en todos los casos evaluados hasta ahora, para proceder a una resolución manual del SEL, en primer lugarse debe representarse la matriz aumentada [A|b]:

$$\begin{align} \begin{bmatrix} a_{11} & a_{12} & . & . & . & a_{1n} &|& b_{11}\\ a_{21} & a_{22} & . & . & . & a_{2n} &|& b_{21} \\ a_{31} & a_{32} & . & . & . & a_{3n} &|& .\\ . & . & . & . & . & . &|& .\\ . & . & . & . & . & . &|& .\\ . & . & . & . & . & . &|& .\\ a_{m1} & a_{m2} & . & . & . & a_{mn} &|& b_{m1} \end{bmatrix} \end{align}$$

Constituida por la matriz de relaciones insumos/recursos por producto/agente (por sitio/situación temporal/condición) y el m-vector columna en donde se especifica el stock total disponible para cada insumo/recurso (en pocas palabras, los datos disponibles). Nuevamente, la estrategia consiste en redudcir a [A|b] por operaciones elementales por fila a su forma escalonada [C|d]. Para el caso que el sistema sea compatible y determinado, se sabrá que x=d.

Modelos Dinámicos

Todo modelo de sistema dinámico representa la evolución de una variable o un conjunto de variables en el tiempo. Si esta evolución depende linealmente de la misma variable o de una combinación lineal del resto de las variables, diremos que la dinámica es lineal

  • Si la dinámica se describe a intervalos regulares (equidistantes en el tiempo), a los cuales se los denomina paso de cálculo, diremos que es un modelo de sistema dinámico en tiempo discreto

  • En consecuencia, la forma más sencilla estaría caracterizada por un sistema de ecuaciones lineales compuesto por una única ecuación: $$ax(k)=x(k+1)$$

en donde $x(k)$ es el valor de la variable de estado $x$ a inicios del $k$-ésimo intervalo de cálculo (o en el instante $k\Delta t$. $\Delta t$ es la duración del paso) y $x(k+1)$ es el valor a inicios del $(k+1)$-ésimo paso de cálculo.

En Ecología y Cs. Ambientales este tipo de modelos se utilizan para la representación de la dinámica de poblaciones.

1.A Caso Unidimensional

Ejemplo 1. Modelo Dinámico Lineal Unidimensional en Tiempo Discreto

A fin de complementar esta descripción teórica puede consultarse este artículo, el cual forma parte de las fichas y apuntes de esta Comisión e incluye material de evaluación.

Supongamos el crecimiento en el tiempo de una población de bacterias, en un caldo de cultivo (aislado y sin intercambio con la frontera exterior). Asumiendo que los recursos son ilimitados, podría sostenerse la hipótesis que establece que esta población presentará una tasa de reproducción con valor constante $r$ (en individuos/individuo) y una tasa de mortalidad con valor constante $m$ (proporción de individuos que sobreviven durante el intervalo de cálculo, computada sobre la base de la población a inicios de este).

Luego, se sabe que el crecimiento de la población será igual a la diferencia entre la cantidad nueva de individuos, por reproducción de las bacterias, y la cantidad de muertes. Si el proceso:

  • Para un muestreo a intervalos regulares de duración $\Delta t=t_2-t_1$ (siendo $t_1$ y $t_2$ dos instancias de medición/observación adyacentes),

  • Si denominamos $x(k)$ a la población de bacterias para el $k$-ésimo intervalo de tiempo (o el instante $t=k \Delta t$)

reagrupando las tasas (coeficientes) se obtiene el siguiente modelo dinámico en tiempo discreto:

$$\begin{align*} rx(k)-mx(k)&=x(k+1) \\ (r-m)x(k)&=x(k+1) \end{align*}$$

análoga a:

$$\begin{align*} ax(k)&=x(k+1) \end{align*}$$

con:

$$\begin{align*} a&=(r-m) \end{align*}$$

Esto es, $a$ representa la tasa de crecimiento efectivo, $x(k+1)$ representa la población de bacterias para el $(k+1)$-ésimo intervalo de duración $\Delta t$ (o instante $t=(k+1) \Delta t$ ). Por esto mismo, puede deducirse que dado un valor inicial $x(k=0)$, también expresado como $x(0)$, y conocidos los coeficientes $r$ y $m$, consecuentemente, para el p-ésimo paso de cálculo, la población de bacterias será:

$$a^px(0)=x(p)$$

En efecto, esta es una propiedad fundamental de los sistemas dinámicos lineales en tiempo discreto. Asimismo, se podrá observar:

  • Si el valor $a>1$, luego la población crece continuamente siguiendo una ley exponencial (superpoblación)

  • Si el valor $a=1$, la cantidad de individuos de la población permanecerá estable (valor constante, la cantidad de nuevos individuos es exactamente la misma cantidad de muertes)

  • Finalmente, si $0 < a < 1$ la población decrecerá exponencialmente con asíntota horizontal 0 (extinción)

Fig.6 Ejemplos de sistemas dinámicos lineales en tiempo discreto, para distintos valores del coeficiente de crecimiento a, representando crecimiento exponencial, estabilidad y extinción/recesión exponencial (a=1.1,1 y 0.8, respectivamente)

Corolario: el valor del coeficiente $a$ brinda la información necesaria para estimar qué tipo de evolución presentará el sistema (superpoblación, estabilidad o extinción).

1.B Caso N-dimensional

  • Supongamos que una población que puede dividirse en $n$-distintos intervalos etarios o estadios ontogenéticos (evolutivos)

  • Supongamos que los recursos disponibles para la población, en principio son ilimitados o no poseen efecto variable sobre la dinámica de población*.

  • Además, supongamos que el objetivo perseguido consista en simular la dinámica demográfica o el tránsito de cohortes (e.g obtener un modelo dinámico de pirámides de población)

    Una cohorte se define como conjunto de individuos que comparten un mismo suceso dentro de un cierto período temporal (e.g. las generaciones)

  • En principio, supondremos que:

    • Entre los pasos $k$ y $k+1$, o lo que es lo mismo entre los instantes $k\Delta t$ y $(k+1)\Delta t$, los individuos sólo pueden:

      • Moverse al estrato o estadio posterior
      • Reproducirse
      • Perecer
    • Asimismo, cada cohorte no puede reproducirse durante su estancia en el primer estrato etario o estadio ontogenéticos

    • Finalmente, cuando una cohorte alcanza el último estrato o estadios ontogenéticos perece en su totalidad durante la duración del paso de cálculo.

Este tipo de modelo lineal $n$-dimensional de la dinámica poblaciones, con estas restricciones, se denomina Modelo de Leslie (ver aquí), en memoria de quien lo formulara por primera vez. Asimismo, dado que en las poblaciones sexuadas la cantidad de hembras es semejante a la de los machos y debido a que las hembras son en efecto quienes brindan nuevos individuos, el modelo suele formularse incialmente para la población de hembras.

Nótese que los supuestos imponen restricciones sobre la duración $\Delta t$ mínima posible para el intervalo de muestreo o paso de cálculo. Esto es, si no se admite permanencia en un estrato etario y la amplitud del rango de edades edades de cada estrato es de $T$ años, y si esta es constante para todos los estratos, puede observarse que debe satisfacerse $T$ $\leq$ $\Delta t$ < $2T$. En caso que los estratos no tengan igual amplitud de rango de edades, podrá observase que deberá satisfacerse $T_1$ $\le$ $\Delta t$ < $2T_2$, en donde $T_1$ y $T_2$ serán respectivamente las amplitudes mínimas y máximas de los distintos intervalos etarios. Asimismo, luego se verá en el modelo propuesto por Lefkovich que, al admitir permanencia en las distintas comoponentes del sistema poblacional, esta restricción se elimina (si bien el sistema de ecuaciones es más complejo).

Modelo de Leslie para una población compuesta por N estratos etarios (ejemplo con N=3)

Para mostrar cómo formular una implementación del Modelo de Leslie, postulemos un ejemplo generalizable para un sistema relativamente simple, compuesto por 3 estratos etarios. De acuerdo a los principios establecidos precedentemente, la estructura del sistema quedaría dada como se muestra en la Fig. 7.

Fig.7 Diagrama de flujo de un modelo de Leslie para un sistema poblacional de 3 estratos. Eneste caso $m_i$ representa la tasa de mortalidad de los individuos presentes en el estrato $x_i$ entre los pasos de cálculo $k$ y $k + 1$. Luego, en cada caso transitará una proporción $(1 − m_i)$ de la población del estrato $x_i$ al estrato $x_{i+1}$, entre el paso $k$ y el paso $k + 1$. Por otro lado, $r_i$ representa la tasa de reproducción de la población del estrato $x_i$ entre $k$ y $k + 1$, dando lugar al inicio del tránsito de una nueva cohorte

Como se sabe, desde la óptica de Teoría de Sistemas y sus aplicaciones en Ecología y Cs Ambientales, entenderemos a todo sistema como un conjunto de componentes (en este caso los estratos etarios), caracterizadas cada una por variables de estado (en este caso, la población situada en cada componente a inicios del $k$-ésimo paso de cálculo), las cuales pueden variar de acuerdo a los flujos de materia/energía entre las componentes (en este caso el tránsito de cohortes, caracterizado por tasas de reproducción y mortalidad).

Por otro lado, para poder formalizar nuestro análisis:

  • Denotaremos como $x_{i}(k)$ a la población del $i$-ésimo estrato a inicios del $k$-ésimo paso de cálculo o intervalo temporal.

  • En principio, denotaremos como $r_i$ a la tasa de reproducción del $i$-ésimo estrato, esto es la cantidad de individuos por individuo que aporta la $i$-ésima componente o estrato al primer estrato o componente.

  • En principio, denotaremos como $(1-m_i)$ a la tasa de supervivencia del tránsito de los individuos situados en el $i$-ésimo estrato (hacia el estrato posterior), durante el $k$-ésimo paso de cálculo o intervalo temporal, siendo $m_i$ la tasa de mortalidad del $i$-ésimo estrato o componente.

A partir de la información provista en la Fig. 6 y bajo los principios asumidos en el Modelo de Leslie, se puede formalizar la expresión matemática que caracteriza la dinámica en cada componente o estrato. Por ejemplo, se sabe que la variación de la población del primer estrato dependerá del aporte realizado por las cohortes que están transitando el segundo y tercer estrato (tasas $r_2$ y $r_3$), ya que en el primer estrato la población no se reproduce. Luego, resulta evidente que para este sistema:

$$\begin{align*} 0x_{1}(k)+r_2x_{2}(k)+r_3x_{3}(k)=x_{1}(k+1) \end{align*}$$

Asimismo, se puede continuar y formular la ecuación que caracteriza la dinámica del segundo estrato. Este estrato recibirá a la población superviviente del primer estrato, por tránsito de cohortes. Esto es:

$$\begin{align*} (1-m_{1})x_{1}(k)+0x_{2}(k)+0x_{3}(k)=x_{2}(k+1) \end{align*}$$

Finalmente, ídem razonamiento podemos utilizar para describir la dinámica del tercer estrato:

$$\begin{align*} 0x_{1}(k)+(1-m_{2})x_{2}(k)+0x_{3}(k)=x_{3}(k+1) \end{align*}$$

Estas son ecuaciones lineales. Por lo que se obtiene un SEL de 3x3:

$$\begin{align*} \begin{matrix} 0x_{1}(k)&+&r_2x_{2}(k)&+&r_3x_{3}(k)&=&x_{1}(k+1) \\ (1-m_{1})x_{1}(k)&+&0x_{2}(k)&+&0x_{3}(k)&=&x_{2}(k+1) \\ 0x_{1}(k)&+&(1-m_{2})x_{2}(k)&+&0x_{3}(k)&=&x_{3}(k+1) \end{matrix} \end{align*}$$

Cuya expresión como producto matricial es:

$$\begin{align*} \begin{bmatrix} 0 & r_2 & r_3 \\ (1-m_{1}) & 0 & 0 \\ 0 & (1-m_{2}) & 0 \end{bmatrix} \begin{bmatrix} x_1(k) \\ x_2(k) \\ x_3(k) \end{bmatrix}= \begin{bmatrix} x_1(k+1) \\ x_2(k+1) \\ x_3(k+1) \end{bmatrix} \end{align*}$$

para simplificar la notación y poder extenderla a casos maś generales, se postula la siguiente igualdad:

$$\begin{align*} \begin{bmatrix} 0 & r_2 & r_3 \\ (1-m_{1}) & 0 & 0 \\ 0 & (1-m_{2}) & 0 \end{bmatrix}= \begin{bmatrix} 0 & a_{12} & a_{13} \\ a_{21} & 0 & 0 \\ 0 & a_{32} & 0 \end{bmatrix} \end{align*}$$

La Matriz de la derecha se denotará como A. Luego, sus coeficientes $a_{ij}$ indicarán el aporte de individuos desde el $j$-ésimo estrato o estadio (columna) que contibuye a la población del i-ésimo estrato o estadio (fila), en cada paso de cálculo o intervalo temporal. Luego, la matriz A es una matriz de transición pues indica la dinámica del tránsito de individuos desde un estrato a otro, representados por las respectivas filas y columnas.

Luego, sea $\mathbf{x}(k)$ el vector columna compuesto por la cantidad de individuos en cada estrato o estadio en el $k$-ésimo intervalo (n-vector de estados del $k$-ésimo intervalo, considerando la población en cada uno de los $n$-estratos o estadios) y sea $\mathbf{x}(k+1)$ el vector de estados para el $(k+1)$-intervalo, este tipo de sistema dinámico queda descrito, en notación simplificada, mediante el producto matricial:

$$\mathbf{Ax}(k)=\mathbf{x}(k+1)$$

Particularmente, en esta implementación el tipo de matriz de transición obtenida recibe el nombre de Matriz de Leslie.

A partir de este modelo $n$-dimenisonal y en comparación con el modelo Unidimensional, puede observarse:

  • La ecuación matricial es semejante a la expresión escalar para el caso unidimensional. En todo caso, A no es una magnitud escalar y no podremos realizar una inferencia direacta como se puede realizar a partir de la tasa a.

  • Si la matriz posee $n$ filas (y por tanto, $n$ columnas), el modelo consta de $n$ estratos o estadios ontogenéticos.

  • En la primera fila se resumen las tasas de reproducción. Si se consideran las filas 2 en adelante, se observa el desarrollo de una diagonal (aquellas entradas con subíndices con $i+1=j$), la cual resume las tasas de supervivencia en el tránsito de una cohorte en el sistema.

  • Finalmente, es posible demostrar que dado un vector inicial de población $\mathbf{x}(0)$, debe cumplirse:

$$\mathbf{A}^{n} \mathbf{x}(0)=\mathbf{x}(n)$$

ya que los coeficientes tienen valor constante.

Finalmente, extendiendo la generalización para cualquier Modelo de Leslie para una población compuesta con $n$ estratos o estadios, tendrá la siguiente estructura:

$$\begin{align*} \begin{bmatrix} 0 & a_{12} & . & . &. & a_{1n} \\ a_{21} & 0 & . &. & . & 0 \\ 0 & . & . &. & . & 0\\ 0 & . & . &. & . & 0\\ 0 & . & . &. & . & 0\\ 0 & 0 & . &. & a_{n(n-1)} & 0 \\ \end{bmatrix} \begin{bmatrix} x_1(k) \\ x_2(k) \\ . \\ . \\ . \\ x_n(k) \end{bmatrix}= \begin{bmatrix} x_1(k+1) \\ x_2(k+1) \\ . \\ . \\ . \\ x_n(k+1) \\ \end{bmatrix} \end{align*}$$

La operación de este tipo de modelos constituye en la aplicación de una transformación matricial lineal $f(\mathbf{x}(k))$ de $R^n\to R^n$ (ver sección 1.5 Transformaciones Matriciales, en Álgebra Lineal de B. Kolman y D. Hill, para mayor información). Específicamente, la matriz de transformación está constituída por la matriz de transición A y dado un vector de estados $\mathbf{x}(k)$, para el $k$-ésimo intervalo, se puede detrminar el valor del vector de estados $\mathbf{x}(k+1)$, mediante el producto matricial entre A y $\mathbf{x}(k)$. En síntesis, la operación del sistema describe la evolución dinámica del vector de estados $\mathbf{x}(k)$, a partir de un vector inicial $\mathbf{x}(0)$, en un sistema de referencias rectanglar (espacio) $n$-dimensional (considerando los $n$ estratos como los ejes).

Modelo de Lekovich para una población compuesta por N estratos etarios (ejemplo con N=3)

El Modelo de Leslie, hemos visto, se formuló inicialmente para poblaciones con estructura de edades (esto es, cada componente está constituída por un intervalo etario, si bien a fin de generalizar se introdujo el concepto de estadio ontológico como componente). Sin embargo, la variable edad no necesariamente es la más importante si el objetivo perseguido es estudiar la evolución de muchos organismos. Por ejemplo, en el caso de los insectos los individuos pasan por distintas etapas o estadios ontogenéticos:

  • Huevo
  • Larva
  • Crisálida / Pupa
  • Adulto

Así, por ejemplo, la tasa de supervivencia puede estar más influenciada por el estadio ontogenético que por la edad en sí misma, lo mismo que la tasa de reproducción (cada etapa está caracterizada por una tasa de supervivencia y de reproducción particulares). Por ejemplo, la supervivencia de un escarabajo no depende de que tenga 3 o 6 meses, sino de que sea una larva o que se encuentre en la etapa adulta.

En 1965, L.P Lefkovich en su artículo The study of population growth in organisms grouped by stages propone un modelo que generaliza la matriz de Leslie, para poblaciones con estructura de tamaño, de estadios, etapas o mixtas. Así, el modelo de Lefkovich admite permanencias en un estadio entre pasos de cálculo o intervalos de tiempo (e.g. la residencia de una cohorte en un estadio puede ser mayor a la duración del $\Delta t$ del $k$-ésimo paso o intervalo de tiempo, los individuos de la primer componentes pueden reproducirse). Luego considerando una estructura de una población compuesta por 3 estadios ontológicos (e.g. como se hizo con los estratos etarios para el caso de Leslie), el esquema de este sistema más generalizable queda dado por el diagrama de flujo que se muestra en la Fig. 8.

Fig.8 Diagrama de flujo de un modelo de Lefkovich para un sistema poblacional de 3 estadios ontológicos. En este caso, se permite que parte de la cohorte permaneza durante todo el intervalo de cálculo dentro del estrato en que se encontraba a inicios del mismo, indicándose mediante las tasas $a_{ij}$, con $j=i$ (la diagonal principal). En otras palabras, se admite un prorateo entre permanencia y movilidad para los supervivientes de la cohorte durante el intervalo de cálculo, dando lugar a considerar que no todos los individuos de la cohorte más anciana perezacan en el mismo. Asimismo, se admite la posibilidad o interpretación que los individuos de la cohorte más joven puedan reproducirse.

Procediendo de forma semejante al caso anterior, se obtendría el siguiente SEL de 3x3:

$$\begin{align*} \begin{matrix} a_{11}x_{1}(k)&+&a_{12}x_{2}(k)&+&a_{13}x_{3}(k)&=&x_{1}(k+1) \\ a_{21}x_{1}(k)&+&a_{22}x_{2}(k)&+&0x_{3}(k)&=&x_{2}(k+1) \\ 0x_{1}(k)&+&a_{32}x_{2}(k)&+&a_{33}x_{3}(k)&=&x_{3}(k+1) \end{matrix} \end{align*}$$

Los valores 0, en las componentes o estadios 2,...,$n$ se corresponden con el hecho que el tránsito siempre es hacia el estadio posterior.

Luego, la expresión del SEL como producto matricial es (verifique):

$$\begin{align*} \begin{bmatrix} a_{11} & a_{12} & a_{13} \\ a_{21} & a_{22} & 0 \\ 0 & a_{32} & a_{33} \end{bmatrix} \begin{bmatrix} x_1(k) \\ x_2(k) \\ x_3(k) \end{bmatrix}= \begin{bmatrix} x_1(k+1) \\ x_2(k+1) \\ x_3(k+1) \end{bmatrix} \end{align*}$$

Pudiéndose generalizar para $n$ estadios ontológicos mediante (verifique):

$$\begin{align*} \begin{bmatrix} a_{11} & a_{12} & . & . &. & a_{1n} \\ a_{21} & a_{22} & . &. & . & 0 \\ 0 & . & . &. & . & 0\\ 0 & . & . &. & . & 0\\ 0 & . & . &. & . & 0\\ 0 & 0 & . &. & a_{n(n-1)} & a_{nn} \\ \end{bmatrix} \begin{bmatrix} x_1(k) \\ x_2(k) \\ . \\ . \\ . \\ x_n(k) \end{bmatrix}= \begin{bmatrix} x_1(k+1) \\ x_2(k+1) \\ . \\ . \\ . \\ x_n(k+1) \\ \end{bmatrix} \end{align*}$$

En donde la diagonal principal agrupa a todos los coeficientes de permanencia.

Análisis de la dinámica de sistemas lineales (Ecología de poblaciones)

Hemos visto que si para una población (sistema), dividida en $n$-estratos/clases/estadios, las tasas de reproducción y mortalidad son constantes (coeficientes), y además se suponen nulos o sin efectos los fenómenos de frontera (i.e. migraciones), mediante la representación gráfica del sistema (identificación/visualización), dando cuenta del fenómeno de tránsito de una cohorte (tasas de supervivencia entre clases/componentes) o de la reproducción de estas (e iniciando nuevos tránsitos), es posible obtener un SEL de nxn, el cual puede representarse mediante el producto matricial:

$$\mathbf{Ax(k)}=\mathbf{x(k+1)}$$

En donde:

  • $\mathbf{A}$ es la matriz de transición, aquella que informa la proporción de individuos que la j-ésima clase aportará a la i-ésima clase (entradas $[a_{ij}]$), entre los intervalos temporales $k$ y $k+1$. Ya sea mediante la reproducción y generación de una nueva cohorte o por el tránsito de una cohorte.

  • $\mathbf{x(k)}$ es el $n$-vector (columna) de población, correspondiente al inicio del $k$-ésimo intervalo.

  • $\mathbf{x(k+1)}$ es el $n$-vector (columna) de población, correspondiente al inicio del $k+1$-ésimo intervalo.

Luego, debido a que las tasas son constantes es fácil demostrar (por inducción) que si se conocen los valores de los elementos de $\mathbf{A}$ y un vector inicial de población $\mathbf{x(0)}$, el vector de población correspondiente al p-ésimo paso queda dado mediante:

$$A^p x(0) = x(p) $$

Así, en primer lugar, podríamos observar que conocidos $\mathbf{A}$ y $x(0)$ es posible operar el sistema incrementando el valor de $p$ y evaluar si la tendencia de $x(p)$ es de crecimiento (exponencial), estabilización hacia un valor constante en cada componente, o decrecimiento (exponencial).

Asimismo, se ha visto que para el caso del sistema lineal más simple (unidimensional), $\mathbf{A}$ se reduce a un valor escalar $a$, el cual es la tasa de crecimiento. Luego, evaluando si $a>$1$, si $a=1$ o si $0&lt;$ $a&lt;1$, puede determinarse de forma analítica y robusta la tendencia del sistema.

En consecuencia, podría advertirse primeramente que esto también debiera poder deducirse de la información presente en la matriz $\mathbf{A}$ de dimensiones $nxn$

Consecuentemente, podríamos formular los siguientes interrogantes a partir de aquello que se viene observando:

  • ¿Existirá alguna magnitud escalar $\lambda$, que pueda computarse a partir de los datos e información presentes en la matriz $\mathbf{A}$ y que permita establecer cuál es la tendencia de crecimiento para un sistema en particular?

  • Además, si se representa a la distribución de población en el k-ésimo intervalo mediante el $n$-vector:

$$\mathbf{\rho(k)}= \begin{bmatrix} \dfrac{x_1(k)}{x_t(k)} \ \\ \dfrac{x_2(k)}{x_t(k)} \ \\ . \\ . \\ . \\ \dfrac{x_n(k)}{x_t(k)} \ \\ \end{bmatrix} $$

  • en donde $x_t(k)$ es la población total del sistema para el k-ésimo intervalo definida mediante la operación:

$$\sum_{i=1}^n{x_i(k)}=x_t (k)$$

De ahí que cada elemento del vector representa a la fracción de la población total que le corresponde a cada clase/estrato/estadio en dicho intervalo (i.e. la distribución de población en dicho intervalo).

Surge otro posible interrogante de interés en Ecología:

  • ¿A partir de algún valor $k$ de intervalo temporal es posible observar que $\mathbf{\rho(k)}$ se estabilice, si bien $\mathbf{x(k)}$ puede continuar creciendo, estable o disminuyendo?

para responder a la primer pregunta, debe introducirse el concepto de autovalor $\lambda$ y su técnica de cómputo. Para responder el último interrogante debe introducirse el concepto de autovector.

Concepto de autovalor (autovector)

Sea $$\mathbf{Ax(k)}=\mathbf{x(k+1)}$$

Un modelo matricial que describe la dinámica poblacional para un sistema de $n$ estratos o clases (etarias u ontogenéticas), con tasas de reproducción y supervivencia constantes, así como eventualmente de permanencia. Si este sistema tiende al crecimiento o decrecimiento monótono o a la estabilidad, luego para un valor arbitrario $k$, relativamente elevado, se podría observar que:

$$\mathbf{Ax(k)}=\lambda \mathbf{x(k)}$$

Igualando ambas expresiones se obtiene:

$$\lambda \mathbf{x(k)}=\mathbf{x(k+1)}$$

Concluyendo que el valor $\lambda$ cumplirá exactamente la misma función que el valor $a$ del caso más simple. Esto es, $\lambda$ representa la tasa constante de crecimiento que el sistema alcanzará a partir de algún momento $\Delta t k$, al respecto del inicio de la operación, aspecto propio la linealidad.

Así, debido a la linealidad presente (i.e. coeficientes de transición constantes) es posible afirmar que en estos sistemas $\lim_{k \to M} \mathbf{Ax}(k) = \lambda \mathbf{x}(k)$. Esto es, a medida que $k$ tiende a aproximarse a un númeror relativamente grande $M$, el producto matricial se transforma en un producto escalar por el factor $\lambda$. Más específicamente, el sistema entrará en régimen a partir de este punto, autodeterminado el crecimiento por el autovalor. Asimismo, podemos pensar que si el ciclo de un tránsito de cohorte es de $T$ años, luego debe cumplirse que $M&gt;T$.

La magnitud $\lambda$ se denomina autovalor dominante del sistema dinámico lineal. Para su cálculo se debe recurrir a principios básicos del Álgebra de Matrices. En primer lugar, se sabe que:

$$\mathbf{Ax(k)}=\lambda \mathbf{x(k)}$$

Luego, realizando las operaciones pertinentes, la igualdad es equivalente a:

$$\mathbf{Ax(k)}-\lambda \mathbf{x(k)}=0$$

Aplicando propiedades fundamentales del Álgebra de Matrices, considerando $\mathbf{x(k)}$ como factor común, se puede postular (verifique):

$$\begin{align*} \mathbf{Ax(k)}-\lambda \mathbf{x(k)}&=\mathbf{I_n}[\mathbf{Ax(k)}-\lambda \mathbf{x(k)}] \end{align*}$$

$$\begin{align*} \mathbf{I_n}[\mathbf{Ax(k)}-\lambda \mathbf{x(k)}]&=\mathbf{I_n}\mathbf{Ax(k)}-\mathbf{I_n}\lambda \mathbf{x(k)} \end{align*}$$

$$\begin{align*} \mathbf{I_n}\mathbf{Ax(k)}-\mathbf{I_n}\lambda \mathbf{x(k)}&=\mathbf{Ax(k)}-\mathbf{I_n}\lambda \mathbf{x(k)} \end{align*}$$

Finalmente, reagrupando e igualando con la expresión precedente, puede demostrarse que:

$$(\mathbf{A}-\mathbf{I_n}\lambda)\mathbf{x(k)}=0$$

Esto constituye un sistema homogéneo $\mathbf{Bx}=\mathbf{0}$, con matriz de coeficientes $\mathbf{B}=(\mathbf{A}-\mathbf{I_n}\lambda)$.

Brevemente, por definición un sistema homogéneo es todo aquel para el cual todos los elementos del vector de términos independientes adquieren valor 0. Asimismo, los sistemas homogéneos presentan la siguiente propiedad:

  • Siempre tienen al menos una solución, puesto que si $\mathbf{x(k)}=(0,0,...,0)=\mathbf{0}$ (todos sus elementos valen 0), se satisface la condición establecida en el sistema (verifique!).

  • Esta solución se denomina solución trivial, puesto que es evidente.

Aun así, un vector de población con todos sus elementos iguales a 0 no tiene significado alguno en la dinámica de poblaciones, al menos que aun se sostenga el principio de generación espontánea, ya refutado en 1861 por el experimento de Pasteur (y ciertamente esta generacióne spontánea no tendría mucho de lineal!).

Por otro lado, se sabe que un sistema compatible puede ser determinado o indeterminado. Evidentemente, los sistemas homogéneos son compatibles. Esto último significa que el sistema puede tener más de una solución, ser indeterminado. Y, además, se sabe que el cómputo del determinante de la matriz de coeficientes del sistema permite inferir esto. En pocas palabras, se sabe que si el valor del determinante es distinto de 0, el sistema es determinado. Por otro lado, si el valor es igual a 0, el sistema es indeterminado.

En caso que el sistema sea indeterminado, no tendrá una solución única, si no un conjunto de vectores solución en $R^n$. Todos los vectores que constituyen esta solución son los denominados autovectores. Luego, este es el caso de interés para nuestra aplicación científica. Así, se postula:

$$det(\mathbf{B})=\mathbf{0}$$

Lo cual equivale a:

$$det(\mathbf{A}-\mathbf{I_n}\lambda)=0$$

Expandiendo la notación, esto es:

$$\left| \begin{array}{cccccc} a_{11} - \lambda & a_{12} & . & . & . & a_{1n} \\ a_{21} & a_{22} - \lambda & . & . & . & a_{2n} \\ . & . & . & . & . & . \\ . & . & . & . & . & . \\ . & . & . & . & . & . \\ a_{n1} & a_{n2} & . & . & . & a_{nn} - \lambda \\ \end{array} \right|=0$$

En donde $|\mathbf{A-I_n\lambda}|$ representa la operación de cómputo del determinante sobre la matriz $\mathbf{B}=\mathbf{A}-\mathbf{I_n}\lambda$. Además se puede notar que la diagonal principal está compuesta por la sustracción de $a_{ii}-\lambda$. Luego, para los modelos más simples que no admiten coeficientes de permanencia (e.g. Modelo de Leslie), esta diagonal quedaría compuesta sólo por elementos con valor $-\lambda$.

Aplicando las propiedades del cómputo de determinantes puede verificarse que el desarrollo del cómputo del determinante en el miembro izquierdo de la ecuación conduce a un polinomio de grado n para $\lambda$, denominado polinomio característico. Finalmente, como el miembro derecho es igual a 0, la labor consiste en hallar las raíces del polinomio característico. En consecuencia, el conjunto solución estará conformado por el vector de valores:

$$\mathbf{\lambda}=(\lambda_1,\lambda_2,...,\lambda_n)$$

Puesto que un polinomio de grado $n$ tiene a lo sumo $n$ soluciones reales.

Además, se define como autovalor dominante, al elemento $\lambda_i$ de mayor valor absoluto. Así si, n=3 y:

$$\mathbf{\lambda}=(-2,0,4)$$

el autovalor dominante, para este caso, es $\lambda_3=4$.

El autovalor dominante será el que brinde la información sobre la tendencia del sistema poblacional.

COROLARIO

Si el autovalor dominante $\lambda&gt;1$ luego el sistema crecerá de forma exponencial. Si el autovalor dominante $\lambda=1$, el sistema tiende a estabilziarse en un vector de población constante. Finalmente si $0&lt;$ $\lambda&lt;1$, el sistema se extingue exponencialmente.

Ejemplos Teóricos

Consideremos un sistema población que se ajusta al modelo de Leslie y está compuesto por 3 clases. En su forma matricial, puede expresarse mediante:

$$\begin{bmatrix} 0 & a_{12} & a_{13} \\ a_{21} & 0 & 0 \\ 0 & a_{32} & 0 \\ \end{bmatrix} \begin{bmatrix} x_{1}(k) \\ x_{2}(k) \\ x_{3}(k) \\ \end{bmatrix}= \begin{bmatrix} x_{1}(k+1) \\ x_{2}(k+1) \\ x_{3}(k+1) \\ \end{bmatrix}$$

Luego, se desarrolla la ecuación $det(\mathbf{A}-\mathbf{I_n}\lambda)=0$, expandiendo esto es:

$$\left| \begin{array}{ccc} -\lambda & a_{12} & a_{13} \\ a_{21} & - \lambda & 0 \\ 0 & a_{32} & -\lambda \\ \end{array} \right|=0$$

Aplicando las propiedades del cómputo de determinantes se obtiene el polinomio característico para este tipo de sistema (Leslie, 3x3, verifique):

$$− \lambda^3 + \lambda a_{12} a_{21} + a_{13} a_{21} a_{32} = 0 $$

Asimismo, para el caso de un sistema poblacional de 3 estadios ontogenéticos que se ajuste al modelo de Lefkovich puede observarse que el polinomio característico queda dado por (verifique):

$$ − \lambda^3 + ( a_{11} + a_{22} + a_{33} ) \lambda^2 + ( a_{12} a_{21} − a_{11} a_{22} − a_{11} a_{33} − a_{22} a_{33} ) \lambda + a_{13} a_{21} a_{32} + a_{11} a_{22} a_{33} − a_{12} a_{21} a_{33} = 0 $$

Finalmente, se procede hallando las raíces o conjunto solución del polinomio característico de grado $n$. El autovalor dominante será el que permitirá establecer si la tendencia es de crecimiento exponencial, estabilización o extinción exponencial. En pocas palabras, el autovalor dominante $\lambda$ es la tasa de crecimiento del sistema evaluado.

Consecuentemente, los pasos a seguir consisten en:

  • Identificar la estructura del Sistema Poblacional (tipo Leslie, Tipo Lefkovich)
  • Formular el sistema dinámico $$\mathbf{Ax(k)}=\mathbf{x(k+1)}$$
  • Obtener el polinomio característico mediante el desarrollo del determinante en la ecuación $$det(\mathbf{A}-\mathbf{I_n}\lambda)=0$$
  • Hallar el conjunto solución del polinomio característico y evaluar el valor del autovalor dominante.
⚠️ **GitHub.com Fallback** ⚠️