Simulating Organogenesis in COMSOL: Parameter Optimization for PDE-based models Denis Menshykau 1,2 , Srivathsan Adivarahan 1 , Philipp Germann 1,2 , Lisa Lermuzeaux 1,3 , Dagmar Iber *1,2 1 D-BSSE, ETH Zurich, Switzerland; 2 Swiss Institute of Bioinformatics (SIB), Switzerland; 3 Department of Bioengineering, University of Nice-Sophia Antipolis, France; *Corresponding author: D-BSSE, ETH Zurich, Mattenstrasse 26, 4058, Switzerland, dagmar.iber@bsse.ethz.ch Abstract: Morphogenesis is a tightly regulated process that has been studied for decades. We are developing data-based and image-basd mechanistic models for a range of developmental processes with a view to integrate the available knowledge and to better understand the underlying regulatory logic. In our previous papers on simulating organogenesis with COMSOL (German et al COMSOL Conf Procedings 2011; Menshykau and Iber, COMSOL Conf Proceedings 2012) we discussed methods to efficiently solve such models on static and growing domains. A further challenge in modeling morphogenesis is the parameterization of such models. Here we discuss COMSOL-based methods for parameter optimization. These routines can be used to determine parameter sets, for which the simulations reproduce experimental data and constraints. Such data is often image based, but may also come from classical biochemical or genetic experiments. Keywords: organogenesis, reaction-diffusion, parameter optimization, COMSOL, image-based modelling. 1. Introduction Increasing computational power and advancement of computational methods now allow us to formulate and solve detailed computational models for gene regulatory networks that control complex biological processes in morphogenesis and organogenesis. The reliability of developmental processes suggests that core processes are deterministic, and accordingly deterministic models for pattern formation have been studied for decades in developmental biology. We formulate the models as systems of reaction- diffusion equations of the form: where denotes the velocity field of the domain and R i the reactions, which couple the equations for the different species X i . D i is the diffusion constant and the Nabla operator. The velocity field might be either imposed, or could be made dependent on the local concentration of proteins, which may affect local cell behavior, e.g. the cell division rate or cell adhesion and motility. A wide range of reaction laws can be used. 1 We have previously shown how to implement and to efficiently solve such models in COMSOL on growing and static multilayer domains. 2,3 The computational models can help to integrate the available knowledge, to test the consistency of current models, and to generate new hypotheses. Previously, we correctly predicted novel genetic regulatory interactions in the limb bud, 4 suggested mechanisms for the control of digit formation, 5 bone development, 6 and for branch point selection in the lung 7 and kidney. 8 One major challenge is the parameterization of computational models. 9 Parameters that are required for the formulation of models include, but are not limited to protein production and degradation rates, diffusion constants, and rate constants for protein-protein interactions. These parameter values can, in principle, be measured in vitro and sometimes also in vivo. However, this is rarely the case and the establishment of measurement methods for a certain system can take years and can typically not be easily transferred to another system. Therefore, direct measurement of all required parameter values is not feasible. Experiments in developmental biology typically deliver information about the organ shape and protein expression patterns at various stages in a WT and mutant system. Here we explore the X i . + ∇(u⋅ X i ) = D i ∇ 2 X i + R i u ∇