Lecture Notes ni lortnoC dna noitamrofnI Sciences Edited yb A.V. Balakrishnan dna .M Thoma 81 I I i I I R gnilledoM dna noitazimitpO fo xelpmoC metsyS Proceed{ngs of eht CT-PIFI ? gnikroW ecnerefnoC Novosibirsk, USSR, 3-9 ,yluJ 8791 Edited yb .G .I Marchuk li I inn I | i • I galreV-regnirpS Berlin Heidelberg New York 1979 Series Editors A. V. Balakrishnan • M. Thoma Advisory Board L. D. Davisson - A. G..1. MacFarlane • H. Kwakernaak Ya, Z. Tsypkin • A. .J Viterbi Editor Prof. Guri Ivanovich Marchuk Computing Center, Novosibirsk 630090, USSR ISBN 3-540-09612-4 Springer-Verlag Bedin Heidelberg NewYork ISBN 0-387-09612-4 Springer-Verlag NewYork Heidelberg Berlin This work is subject to copyright. All the whether whole rights are reserved, or part of the material is concerned, specifically those of translation, -er printing, of re-use illustrations, broadcasting, reproduction photocopying by machine or similar means, dna storage in data .sknab § Under 54 of the Copyright German where Law copies are fomra de other than private use, fee a is payable to the the amount publisher, of the fee to with the agreement determined be by publisher. © gadeV-regnirpS Berlin Heidelberg 979(cid:127) Printed in Germany Printing dna binding: Beltz Offsetdruck, Hemsbach/Bergstr. 2060/3020-543210 PREPACE These Proceedings contain most of the papers presented at the IFIP TC-7 Working Conference on Modelling and Optimization of Complex Systems held in Novosibirsk on 3-9 July, 1978. The Conference was organized by the IPIP Technical Committee on Optimization with the Computing Center of the Siberian Branch of the USSR Academy of Sciences and sponsored by the USSR Academy of Scien- ces. The Conference was attended by 70 scientists from 10 countries. The program offered a broad view of optimization techniques currently in use and under investigation. Major emphasis was on recent advances in optimal control and mathematical programming and their application to modelling, identification and control of large systems, in parti- cular, recent applications in areas such as biological, environmental and sccio-economio systems. The International Programme Committee of the Conference consisted of: A.V.Balakrishnan (Chairman, USA), C.Bruni (Italy), J.L.Lions (France), K.Malanowskl (Poland), G.I.Marchuk (USSR), R.R.Mohler (USA), L.S.Pontryagin (USSR), J.Stoer (FRG). VJ TABLE OF CONTENTS Invited Speakers A.BERTUZZI, C.BRUNI, A.GANDOLFI, G.KOCH ~aximum Likelihood Identification of an Immune Response Model ... I M. JILEK, P.KLEIN Stochastic Model of the Immune Response ......................... 15 A. KALLI AUER A Computational Method for Hierarchical Optimization of Complex Systems . . . . ................................................... . 26 K. D. MALANOWSKI On Regulerity of Solutions to Constrained Optimal Control Problems for Systems with Control Appearing Linearly ........... @I R. R. MOHLER Bilinear Control Structures in Immunology ....................... 58 ~ATHE~ATICAL MODELS IN IMMUNOLOGY A. L. ASACHENKOV, G. I. N~RCHUK Mathematical Model of a Disease and Some Results of Numerical Exp eriment s 69 @ * e * e ~ . * e • . • o • • @ e o • . o . t m e oee.@eo4Be6bm • • ~ Omo • @e ovoes@ L.N. BEL~H, G. I KUHCRAM. Chronic Forms of a Disease and Their Treatment according to Mathematical Immune Response ~odels . 79 o o o o o a o e o v o o e o o 6 1 e 4 J o e o o o o o o B.F.DIRBOV, M.A.LIVSHITS, M.Y.VOLKENSTEIN The Effect of a Time Lag in the Immune Reaction ................ 87 V E. V. GRUN TENKO Immunological Aspects of Neoplastic Growth .................... 95 V. P. LOZOVOY, S .~. SHERGIN Structure-Functional Arrangement of the Immune System ........ 103 I. G. NARCHUK Mathematical ~odels in Immunology and Their Interpretations .. 114 V. R. PETROV Clinical Immunology: Problems and Prospects ................... 130 OPTIMIZATION OF COMPLEX SYSTEMS V. M. ALEXANDROV Approximate Solution of Optimal Control Problems ............. 147 K. BELLMANN, J. BORN Numerical Solution of Adaptation Problems by means of an Evolution Strategy ........................................... 157 V.L.BERESNEV, E.KH.GIMADI, V.T.DEMENTIEV Some Models of Taking Optimal Decisions in Standardization ... 168 J. DOLEZAL On Optimal Control of Discrete Systems with Delays ........... 181 YU. P.DROBYSHEV, V.V.PUKHOV Analysis of the Influence of a System on Objects as a Problem of Transformation of Data Tables ............................. 187 E.V.ILJIN, A.V.MEDVEDEV, N.P.NOVIKOV On Non-Parametric Algorithms of Optimization ................. 198 M.LUCERTINI, A.MARCHETTI SPACCAMELA Approximate Solutions of an Integer Linear Progrsmm~ng Problem with Resource Variations ..................................... 207 IV V.M.NATROSOV, S.N.VASSILYEV. 0.G.DIVAKOV, A.I. TYATUSHKIN On Technology of Modelling and Optimization of Complex Systems • 220 G. I. NARCHUK, V.V. PENENK0 Application of Optimization Methods to the Problem of Mathema- tical Simulation of Atmospheric Processes and Environment ..... 240 A. NARROCC0, 0.PIRONNEAU Optimum Design of a Magnet with Finite Elements .............. 253 A. MASLOWSKI Finite Elements Approximation in State Identification of Large Scale Linear Space-Distributed Systems ........................ 264 A.A.PETROV, I.G.POSPELOV ~athsmatical Modelling of Socio-Economic System .............. 274 M. V. YAKOVLEV On one Representation of Necessary Conditions of Optimality for Discrete Controlled Systems .................................. 284 MAXIMUM LIKELIHOOD IDENTIFICATION OF AN IMMUNE RESPONSE MODEL )$( Alessandro Bertuzzi(1),Carlo Bruni (2'I) Alberto Gandolfi )I( , Giorgio Koch (3,1) )I( Centro di Studio dei Sistemi di Controllo e Calcolo Automatici del CNR, Via Eudossiana 18, 00184 ROMA, ITALY )2( Istituto di Automatica dell'Universit~ di Roma, Via Eudossiana 18, 00184 ROMA, ITALY )3( Istituto Matematico "G.Castelnuovo" dell'Universita di Roma, Citt~ Universitaria, Piazzale delle Scienze 5, 00185 ROMA,ITALY. I. INTRODUCTION In a previous work a model for the humoral immune response was established, based on the clonal selection theory (Bruni et a1.1974, 1975). This model was tested against available data in the literature (Eisen and Siskind, 1964; Siskind et al.,1968;Werblin et al., 1973) giving the parameters some reasonable values, with rather satisfactory results. In order to achieve a more thorough validation of the model, it appeared necessary to a) obtain a more complete and homogeneous set of experimental data, b) achieve a better knowledge of the values of the model parameters. As far as the first point is concerned, an ad hoc experimental program was carried out, in which the population of antibodies gene- rated in the primary and secondary immune response triggered by ino- culation of a suitable antigen was tested at different times (Oratore, 1977). Inbred Guinea pigs (strain 13) were inoculated by doses of 0.1 mg and I mg of DNP-RNAase in Freund's complete adjuvant. Bleedings were performed weekly and the adopted titration procedure was the Farr technique.Also, some experiments were made to get information about the time evolution of the free antigen concentration in the blood during the maturation of the immune response (Oratore, unpubli- shed results). The aim of this paper is to deal with the second point, namely to identify the unknown parameters of the model. )*( This work was partially supported by C.N.R.. The experiments which are referred to in this paper were conducted at Oregon Regional Primate Research Center, under a joint scientific program between the Centro di Studio dei Sistemi di Controllo e Calcolo Automatici and the Oregon State University, supported by Italian National Re- search Council and the U.S. National Science Foundation ( Grant n. ENG-74-15330 AOI). 2. AN IDENTIFICATION ORIENTED MODEL A first problem was to revise the previously proposed model in order to get an identification oriented version of it. This prima- ri~y led to assume as input and output variables the free antigen concentration and respectively the antibody concentration density,i.e. those quantities which were experimentally observed. With reference to the original model, this means to consider only the phenomena of stimulation, differentiation and duplication of lymphocytes, and of antibody production and removal, while the feedback phenomena of antigen removal by the antibodies was cut off (see Fig. 1). In addi- tion, the term which in the original model described a burst of anti- body synthesis by immunocompetent cells immediately afteerh eir duplica tion was disregarded, since it was found to be not quantitatively relevant. By this way we arrived at the following simplified model equations: (K,t) &C =ec 1 KH(t) a-----~ I+KH(t) PS (KH(t))C(K't)-~c I+KH(t) Ps (KH(t))C(K't) + (2.1) 1 C(K,t) + ~Pc(K) c T ~Cp ,K( )t KH )t( _ I ~t =2~c I+KH(t) Ps (KH(t))C(K't) ~p Cp(K,t) (2.2) ~ST(K't) _1 KH(t) 1 I ~(K,t) (2.3) ~t ~sCp(K't)+~s C(K't) B T I+KH(t)ST (K't) s r I+KH(t) where: K is the association constant between antigen and antibody sites ranging in the interval [0,~) (antigen is considered functionally univalent) H(t) is the free antigen concentration in the circulating fluids C(K,t) is the concentration density, with respect to K, of immuno- competent cells Cp(K,t) is the concentration density of plasmacell~ ST(K,t) is the total (free and bound) antibody sites concentration density Ps(KH(t)) is the probability that a cell of affinity K is stimulated; in (Bruni et al., 1975) it is shown that this is a smooth window function, which takes value very ~loseto I if I-~, !s --< KH(t) ~I~ and very close to zero otherwise Pc(K) is the original distribution of immunocompetent cells and it is assumed to be lognormal and known. The parameters which appear in the model, are: c a proliferation rate constant of stimulated immunocompetent cells 8 rate of production of new immunocompetent cells from stem cells Tc mean lifetime of immunocompetent cells T mean lifetime of plasmacells P ~B mean lifetime of the immune complex s T mean lifetime of free antibody sites s a rate constant of antibody production by plasmacells s a ! rate constant of antibody production by immunocompetent cells oi,o2 endpoints of the interval of occupied receptor site ratio causing stimulation. To these parameters we must add suitable initial conditions for equa- tions (2.1)-(2.3). ~ile the input H is directly measured at different times, the total antibody sites concentration density is indirectly measured through the concentration Y(X,t) of bound antibody sites in the titr~ tion assay made on the circulating fluids at different times and for fixed free hapten concentration X: [=. XK S T(K,t) = Y(X,t) j I+KX I+KH(t) dK (2.4) 0 Noting that the free antibody site concentration density in the cir- ST(K,t) culating fluids is given by the term I+KH(t) , the previous relation provides under mild assumptions (Bruni et al., 1976) a one to one correspondence between the free antibody site concentration density itself at time t and the behaviour of Y(X,t) as a function of X. Thus Y(X,t) plays the role of actual model output (Fig. I). prolifera-ll~ TS (K,t) , Stimulation and I~li(t) IAntigen tion of lymphocytes, pro-ll I~=K-'~) I Titration ! I Y(X,t) -I of 1=I removal I duction and removal assay I I antibodies 1=I T .I I Fig. I - Block diagram showing connection among quantitSes relevant for identification purpose. The dashed square includes the identification oriented immune response model; ~i(t) denotes the rate of injection of antigen in the organism. NOW the second problem is to choose those parameters in the model which are to be identified. Indeed, if all parameters and initial conditions appearing in the model were assumed to be unknown in the identification procedure, then we would face: )a difficulties in meeting possible a priori (or experimental) information about rela- tionships among some of them, b) possible non identifiability problems due to low sensitivity of the output with respect to some parameters and/or equivalent effect on the output itself, )c an excessive comput~ tional burden. With the purpose of overcoming these difficulties, only the primary response was considered, so that the initial condi- tions reduce to: C(K,0) = CoPc(K) = rc~Pc(K) Cp(K,0) = 0 (2.5) ST (K,0) = TSe s 'C (K,0) In the first equation (2.5), o C denotes the initial total number per unit volume of immunocompetent cells which turns out to be equal to ~c ,8 since C(K,0) must be a stationary solution of (2.1) in ab- sence of stimulation.By simulation of the model it appeared that variations in o C may be (output) compensated by suitable variations of s a , a'.s Therefore C0 was fixed to the biologically reasonable value: C = 8-10 -16 moles/liter O Information reported in the literature (Weigle, 1961; Gowans, 1970; Mattioli and Tomasi, 1973) funther suggest to assume the following values: = 6 days = 144 h P ~c = 300 days = 7200 h s = B T 2 T B = Co/~ c =1.1.10 -19 moles/liter.h Finally, the parameters which are considered as unknown in the identi- fication procedure are : s'TB This selection seems to meet some remarks (Mohler, 1978) about the sensitivity of the model with respect to various parameters.