A lack of detailed knowledge of the gene regulatory networks governing differentiation of embryonic stem cells (ESC) significantly hinders the development of robust experimental protocols. An integrated experimental and computational approach will have significant impact in developing mature, functional cell types thus advancing the fields of stem cell therapy and regenerative medicine. We have developed a novel algorithm to capture the dynamics of gene regulatory networks involved in the differentiation of ESCs towards the pancreatic lineage from a restricted number of input data set. Our methodology can also capture the effect of external perturbations on the regulatory network. The algorithm employs a bi-level optimization method where the upper level is formulated as an integer programming problem optimizing the network topology, while the lower level deals with continuous variables optimizing the strength of the network connectivity. The input data is the time profile of gene expression during the pancreatic differentiation process. The pancreatic lineage is induced by first differentiating the ES cells to endoderms, followed by subsequent manipulations of extracellular matrix and chemical factors. The algorithm is further validated experimentally by silencing Foxa2 gene. The predicted network has an excellent agreement with the experimental data, as also with the current understanding of pancreas development. This presented integrated methodology holds great promise in developing robust in-silico differentiation protocols.