bioprocess_model.py

# bioprocess_model.py

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from scipy.optimize import curve_fit
from sklearn.metrics import mean_squared_error
from sympy import symbols, lambdify, sympify, Function

class BioprocessModel:
    def __init__(self):
        self.params = {}
        self.r2 = {}
        self.rmse = {}
        self.models = {}  # Initialize the models dictionary

    @staticmethod
    def logistic(time, xo, xm, um):
        return (xo * np.exp(um * time)) / (1 - (xo / xm) * (1 - np.exp(um * time)))

    @staticmethod
    def substrate(time, so, p, q, xo, xm, um):
        return so - (p * xo * ((np.exp(um * time)) / (1 - (xo / xm) * (1 - np.exp(um * time))) - 1)) - \
               (q * (xm / um) * np.log(1 - (xo / xm) * (1 - np.exp(um * time))))

    @staticmethod
    def product(time, po, alpha, beta, xo, xm, um):
        return po + (alpha * xo * ((np.exp(um * time) / (1 - (xo / xm) * (1 - np.exp(um * time))) - 1))) + \
               (beta * (xm / um) * np.log(1 - (xo / xm) * (1 - np.exp(um * time))))

    def set_model(self, model_type, equation, params_str):
        """
        Configures the model based on the type, equation, and parameters.

        :param model_type: Type of the model ('biomass', 'substrate', 'product')
        :param equation: The equation as a string
        :param params_str: Comma-separated string of parameter names
        """
        t_symbol = symbols('t')
        X = Function('X')  # Definir 'X(t)' como una función simbólica

        try:
            expr = sympify(equation)
        except Exception as e:
            raise ValueError(f"Error al parsear la ecuación '{equation}': {e}")

        params = [param.strip() for param in params_str.split(',')]
        params_symbols = symbols(params)

        # Extraer símbolos utilizados en la expresión
        used_symbols = expr.free_symbols
        # Convertir símbolos a strings
        used_params = [str(s) for s in used_symbols if s != t_symbol]

        # Verificar que todos los parámetros en params_str estén usados en la ecuación
        for param in params:
            if param not in used_params:
                raise ValueError(f"El parámetro '{param}' no se usa en la ecuación '{equation}'.")

        if model_type == 'biomass':
            # Biomasa como función de tiempo y parámetros
            func_expr = expr
            func = lambdify((t_symbol, *params_symbols), func_expr, 'numpy')
            self.models['biomass'] = {
                'function': func,
                'params': params
            }
        elif model_type in ['substrate', 'product']:
            # Estos modelos dependen de biomasa, que ya debería estar establecida
            if 'biomass' not in self.models:
                raise ValueError("Biomasa debe estar configurada antes de Sustrato o Producto.")
            biomass_func = self.models['biomass']['function']
            # Reemplazar 'X(t)' por la función de biomasa
            func_expr = expr.subs('X(t)', biomass_func)
            func = lambdify((t_symbol, *params_symbols), func_expr, 'numpy')
            self.models[model_type] = {
                'function': func,
                'params': params
            }
        else:
            raise ValueError(f"Tipo de modelo no soportado: {model_type}")

    def fit_model(self, model_type, time, data, bounds=([-np.inf], [np.inf])):
        """
        Fits the model to the data.

        :param model_type: Type of the model ('biomass', 'substrate', 'product')
        :param time: Time data
        :param data: Observed data to fit
        :param bounds: Bounds for the parameters
        :return: Predicted data from the model
        """
        if model_type not in self.models:
            raise ValueError(f"Model type '{model_type}' is not set. Please use set_model first.")

        func = self.models[model_type]['function']
        params = self.models[model_type]['params']

        # Definir la función de ajuste
        def fit_func(t, *args):
            try:
                y = func(t, *args)
                return y
            except Exception as e:
                raise RuntimeError(f"Error en fit_func: {e}")

        # Definir una estimación inicial para los parámetros
        p0 = [1.0] * len(params)  # Puedes ajustar estos valores según sea necesario

        try:
            # Definir los límites correctamente
            lower_bounds, upper_bounds = bounds

            # Ajustar el modelo usando curve_fit con p0
            popt, _ = curve_fit(fit_func, time, data, p0=p0, bounds=(lower_bounds, upper_bounds), maxfev=10000)

            # Guardar los parámetros ajustados en el modelo
            self.params[model_type] = {param: val for param, val in zip(params, popt)}
            y_pred = fit_func(time, *popt)
            self.r2[model_type] = 1 - (np.sum((data - y_pred) ** 2) / np.sum((data - np.mean(data)) ** 2))
            self.rmse[model_type] = np.sqrt(mean_squared_error(data, y_pred))
            return y_pred
        except Exception as e:
            raise RuntimeError(f"Error while fitting {model_type} model: {str(e)}")