Lecture світів-07

Creation of Worlds — 07. Parameter Selection and Solution Search

Spreadsheets provide a useful capability: they allow for the approximation of relationships describing empirical data, as well as, in the case of MS Excel, – finding the best solutions for describing systems with many degrees of freedom.

7. Parameter Selection and Solution Search 7.1. Parameter Selection A powerful tool in spreadsheets is the parameter selection feature. We will explain its logic first with a simple example, and then show its application to solve a problem where its potential is better demonstrated. The standard parameter selection tool in LO Calc ("Tools / Goal Seek...") allows you to solve the following problem: determine at what value of one cell (variable) the value of another cell (target) will be equal to a certain value. Naturally, this is possible if the target cell contains a formula that depends on the value of the variable cell. Consider the following example. We have a series of numbers for which we know (or have decided) what dependency describes them, but we do not know the value of the parameter included in this dependency. We want to perform an approximation (search for an approximation, an optimal description) of this dependency. Let's use, for example, the exponential growth model, the very first model from our course (sections 2.1 – 2.4). Using our model, we will create a series of numbers describing exponential growth, and then introduce stochastic changes into the obtained result (e.g., within 10% of the value itself in both directions). To do this, it is enough to multiply each value by the expression (1-(0.5-RAND())/5); this expression produces a value uniformly distributed between 0.9 and 1.1. Perhaps you are not sure that the formula (1-(0.5-RAND())/5) produces a distribution between 0.9 and 1.1? You can check this in two or three minutes (if you know how to use spreadsheets). Let's look at this example in more detail, as it is a good opportunity to get acquainted with new modeling tools. Let's create a helper file. Create 1000 cells in it, containing the formula (1-(0.5-RAND())/5). Create another 1000 cells where we calculate the rounded value of this random variable to two decimal places. You understand that it is enough to enter these formulas into two cells and then "drag" them to the remaining 999 pairs of cells? To see such formulas, you can use the method we explained when discussing the binomial distribution in section 6.2. On the left side of Fig. 7.1, formula viewing is enabled (using the "View / Show Formula" command), on the right, the results of their work are demonstrated. Fig. 7.1. Example of checking the distribution produced by a certain formula. This figure combines two ways of viewing models; on the left – viewing formulas ("View / Show Formula"), on the right – viewing values (and a histogram). We can check the minimum and maximum values produced by the formula we are using with the MIN(range) and MAX(range) functions; the use of such formulas is quite intuitive. To assess the nature of the distribution, we should calculate the frequencies of individual classes; we will use the COUNTIF function for this. Let's build a series of numbers from the minimum to the maximum value with a step corresponding to the rounding step; it is clear that it is enough to enter the first two values and "drag" it to the desired range. Next, we will "drag" the cells with the counting function (of course, with absolute addressing where needed). We build the distribution... and see that the deviations from uniformity are most likely random. Do you understand why the values 0.9 and 1.1 should occur, on average, half as often as all other values? Only values "from above" (from 0.9 to 0.905) are rounded to 0.9, only "from below" (from 1.095 to 1.1) are rounded to 1.1, and all other values are rounded both "from above" and "from below" (e.g., 1 is rounded from 0.995 to 1.005). Let's return to the main task of this section: using the formula we just checked for the random variation of the terms in the exponential growth series. In Fig. 7.2, you can see that this random correction can, for example, lead to the fact that the next term in this series is less than the previous one (note the 15th and 16th terms of the series, on lines 21 and 22). Fig. 7.2. Column C shows the values from column B, modified by multiplication by a random variable within ± 10% (note the row for formulas). Assume that based on the values in column C, we need to determine the parameters by which the sequence in column B was constructed (and the values in column C themselves can be considered measurements of the values in column B with a certain error). We need to construct another sequence, which we will build according to our chosen pattern and "fit" (modify it to bring it as close as possible) this pattern to the values in column C. We will do this in column E. Let's fill column E with the same sequence as column B, only instead of the value r, which is unknown to us (for column C), we will enter the value ~ r – an approximation of r (in cell E3). We can set ~ r to anything – for example, to 1, as done in Fig. 7.3. Now we can compare the values we have in column C with the approximations in column E. To do this comparison, in column F, for each row, we will calculate the square of the difference between the corresponding values from C and E. In cell E4, we will calculate the sum of all such squared differences between the values from the two series. Fig. 7.3 shows the formulas used to do this. It can be said that column E contains a model of the original that is in column C; column C is an artificially distorted model of the series in column A. Thus, column E contains an indirect model of the sequence in column A. Fig. 7.3. Construction of a sequence that, through an optimization procedure, should approximate an exponential sequence, where the parameter r is "unknown" for columns other than column B. The right side shows an example of a formula that calculates an exponential sequence based on the value of ~ r; the right side shows an example of a formula that calculates the squared difference between "measurements with error" and our approximation. Note an important circumstance. We are building an approximation to a series of "measurements with error" that follows an exponential law. Because of this, we plan to find the parameter r by which the "unknown" true exponential sequence was constructed. We hope that since there are quite a few terms in this sequence (in Fig. 7.2, you can see that we are considering the first 50 values of the exponential sequence), the errors will balance each other out, and we will find the "true" value of r quite accurately. Is it important for us from which value of ~ r we start our search? No! What is important is not the value from which we start, but how we will approach the quantity we are looking for. By the way, in complex optimization problems with many parameters, the result of the search depends on the initial value, but in this example, we can start with any value (check it yourself!). How to determine how well we have obtained the approximation? Calculate the measure of discrepancy between the model and the original (in our case, exponential growth calculated based on ~ r and the series of "measurements with error"). We will calculate this measure in cell E4 (in Fig. 7.4, you can see that if the value of ~ r is taken from the "ceiling" of the series, the discrepancy measure will be very large). In Fig. 7.4, this discrepancy measure is denoted by a formula, but we can use the notation G (from the English word "gap") for it. In addition, Fig. 7.4 shows how to call the "Tools / Goal Seek..." function. Fig. 7.4. Cell E4 calculates G – the measure of discrepancy between the model and the original; it shows how to call the "Goal Seek..." dialog box. A dialog box will appear in which we will need to specify the target cell with the function that depends on the cell with the variable value. Note that we want to get a value of 0 in the target cell (this is, of course, impossible, as it corresponds to the case where the values in the two columns are identical; in any case, LO Calc will try to approach this value). The variable cell is the cell with the value ~ r (Fig. 7.5). Fig. 7.5. LO Calc will minimize (approach to 0) the value in cell E4 by changing the value in cell E3. Naturally, the value in the target cell cannot reach 0. LO Calc finds the best approximation and reports it (Fig. 7.6). Fig. 7.6. LO Calc suggests accepting the best approximation it found. After the found approximation is accepted, the program will substitute it into the cell for ~ r. It should be noted that the approximation turned out to be quite good (Fig. 7.7). Fig. 7.7. Result of parameter selection. The approximation (0.1411) differs from the value that should have been set (0.14), but the correspondence is quite good. Can this result be improved? Without a doubt. One way is as follows. As the measure G (discrepancy between the model and the original), we used the sum of squared deviations: G 1 =∑((N- ~ N) 2). However, other solutions are possible. The series we are working with is accelerating. The first term of the series is 1, the fiftieth is 700; a 1% deviation from the first term corresponds to 0.01, and from the fiftieth – 7. How to reduce the influence of "scale"? Check for yourself how the accuracy of estimating r will change if you minimize the value G 2 =∑((N- ~ N) 2 /(N+ ~ N))! 7.2. Solving Problems with Multiple Parameters LO Calc also has a tool for finding optimal solutions to problems with multiple variables. This is the "Solver...", a link to which can be found in Fig. 7.4 under the link to "Goal Seek..." in the "Tools" menu. Unfortunately, this LO Calc tool is weaker than the analogous tool in MO Excel, as it only works with linear equations (i.e., it can only be used in relatively simple cases). However, in reality, LO Calc is a more powerful program, as it has an NLPSolver extension, which is suitable for nonlinear problems. To consider the operation of tools for solving problems with multiple variables, we will use data on the age and body length of 100 female lake frogs, Pelophylax ridibundus, taken from the dissertation work of O. Ye. Usova (2016). Table 7.1. Age (in years) and body length (in 0.1 mm) of 100 female lake frogs, Pelophylax ridibundus (Usova, 2016)

Women

Can this result be improved? Undoubtedly. One of the ways is as follows. As the measure G (discrepancy between the model and the original), we applied the sum of squared deviations: G1=∑((N-~N)2). However, other solutions are possible. The series we are working with grows with acceleration. The first term of the series is 1, the fiftieth is 700; a discrepancy of 1% relative to the first term corresponds to 0.01, and relative to the fiftieth corresponds to 7. How to reduce the influence of "scale"? Check independently how the accuracy of the r estimate will change if we minimize the value G2=∑((N-~N)2/(N+~N))!

Women

Can this result be improved? Undoubtedly. One of the ways is as follows. As the measure G (discrepancy between the model and the original), we applied the sum of squared deviations: G1=∑((N-~N)2). However, other solutions are possible. The series we are working with grows with acceleration. The first term of the series is 1, the fiftieth is 700; a discrepancy of 1% relative to the first term corresponds to 0.01, and relative to the fiftieth corresponds to 7. How to reduce the influence of "scale"? Check independently how the accuracy of the r estimate will change if we minimize the value G2=∑((N-~N)2/(N+~N))!

Women

Can this result be improved? Undoubtedly. One of the ways is as follows. As the measure G (discrepancy between the model and the original), we applied the sum of squared deviations: G1=∑((N-~N)2). However, other solutions are possible. The series we are working with grows with acceleration. The first term of the series is 1, the fiftieth is 700; a discrepancy of 1% relative to the first term corresponds to 0.01, and relative to the fiftieth corresponds to 7. How to reduce the influence of "scale"? Check independently how the accuracy of the r estimate will change if we minimize the value G2=∑((N-~N)2/(N+~N))!

Women

Can this result be improved? Undoubtedly. One of the ways is as follows. As the measure G (discrepancy between the model and the original), we applied the sum of squared deviations: G1=∑((N-~N)2). However, other solutions are possible. The series we are working with grows with acceleration. The first term of the series is 1, the fiftieth is 700; a discrepancy of 1% relative to the first term corresponds to 0.01, and relative to the fiftieth corresponds to 7. How to reduce the influence of "scale"? Check independently how the accuracy of the r estimate will change if we minimize the value G2=∑((N-~N)2/(N+~N))!

Women

Can this result be improved? Undoubtedly. One of the ways is as follows. As the measure G (discrepancy between the model and the original), we applied the sum of squared deviations: G1=∑((N-~N)2). However, other solutions are possible. The series we are working with grows with acceleration. The first term of the series is 1, the fiftieth is 700; a discrepancy of 1% relative to the first term corresponds to 0.01, and relative to the fiftieth corresponds to 7. How to reduce the influence of "scale"? Check independently how the accuracy of the r estimate will change if we minimize the value G2=∑((N-~N)2/(N+~N))!

4

245

3

347

5

462

7

644

6

859

3

263

5

349

3

464

5

650

7

861

3

291

3

358

4

486

7

672

5

865

3

304

4

358

5

502

4

685

9

867

3

305

3

359

4

522

5

689

7

868

4

306

4

373

7

524

6

702

6

871

4

310

5

376

5

549

5

715

6

878

3

312

3

378

6

549

9

722

7

883

3

313

4

379

3

557

5

729

7

895

3

319

3

383

6

560

6

743

7

901

5

320

5

384

6

566

7

755

7

903

3

322

5

395

7

568

7

759

10

931

3

324

3

396

6

591

8

785

8

934

3

326

4

408

4

595

8

790

10

959

4

326

3

410

6

596

5

792

8

979

4

335

4

420

7

607

5

801

7

990

4

339

5

424

7

619

7

805

6

1026

4

342

5

425

6

620

8

822

7

1042

5

344

6

434

6

634

8

828

7

1055

5

345

5

447

6

642

6

852

8

1211

Our task is to establish the pattern by which the size of these frogs increases with age. In other words, we will find the regression equation (from Latin regressio – 'backward movement', meaning a transition not from an equation to values, but vice versa – from values to equations). We will use modeling for this. Let's denote the age of individuals by the letter A, and the length by L. Keep in mind that in models describing individual growth, some notations may have a different meaning than in models describing population dynamics (there aren't enough letters in the alphabet!). The simplest mathematical model that can be used to describe individual growth is a straight line: y = a + b*x. The coefficient b determines the angle at which this line passes (for example, if b=1, this line passes at a 45° angle counterclockwise relative to the x-axis). The coefficient a determines the upward or downward shift of this line relative to the origin (if a=1, this line passes through the origin). Let's model the growth of frogs, the data for which are given in Table 7.1, using such an equation. The approach we will use will be analogous to the previous case. We have data on A and L, which depend on A. This is empirical data; it does not depend on us. We will model L using a variable ~ L (L with a tilde), which will depend on coefficients a and b: ~ L = a + b*A. To do this, we need to set specific values for coefficients a and b and calculate for each individual, depending on its age, the expected body length ~ L according to the model. The model will be good if, for all frogs, the empirical L and the calculated ~ L correspond well to each other. To ensure this correspondence, we will need to determine a measure of the differences between the model and the empirical data G, and then minimize G by changing a and b. How to do this is shown in Fig. 7.8. Fig. 7.8. Shows the formulas needed to calculate ~ L and G. To minimize G by changing a and b, you need to install the NLPSolver add-on. How to do this depends on the operating system and other software installed on your computer. For example, on Ubuntu, which the author of this course uses, it is sufficient to execute the following commands: sudo apt update sudo apt install libreoffice-nlpsolver On Windows computers, solving this problem may be more complex (thanks to Kateryna Nesterenko, who figured this out!). To install NLPSolver on such a computer, you need to download the installation file and run it, agreeing to all prompts during installation. If you encounter the error "Could not create Java implementetion loader", or a message about needing a jre environment, you should download the Java installer from here, install it on your computer, and then install NLPSolver. We will assume that you have successfully installed NLPSolver. Call "Solver" (how to do this is shown in Fig. 7.4), set the optimization conditions (Fig. 7.9). Note that you can only solve this problem with NLPSolver installed. Fig. 7.9. Target and variable cells are set. To select an algorithm, go to the "Parameters" menu. In the "Solver Mechanism" dialog, you can choose one of four algorithms (Fig. 7.10). Fig. 7.10. If you have NLPSolver installed, you can choose one of four algorithms. Let's compare these algorithms. Let's propose initial values: a=1; b=1. Both when using "DEPS Evolutionary Algorithm" and when using "SCO Evolutionary Algorithm" with default parameters, the first run of the solver gives a=4; b=4. This is because these algorithms search for solutions at a distance from the initial value that does not exceed a certain limit. Both parameters can change by no more than 4 times. The second run ends with values a=16; b=16, the third with a=64; b=64, and the fourth with a=21,911; b=105,836. The last value does not change with subsequent solver runs. Attempting to select "LibreOffice Linear Solver" does not end well: the program reports that this case is nonlinear. If you select "LibreOffice Swarm Nonlinear Solver (Experimental)", the solver pauses for a while and then immediately finds the solution a=21,911; b=105,836. To understand what we have obtained, let's add a chart. Select the data in columns A, B, and C. Go to "Insert / Chart... / XY Chart / Lines and Points". Select the options "Data rows in columns" and "First column as labels". Add the necessary labels. We will get something similar to the chart shown in Fig. 7.11. Fig. 7.11. To make this chart understandable, you should leave only the icons for the empirical data and only the line for the calculated dependence. After some adjustments (enter the chart, right-click on the line, select "Format Data Series..." and edit as needed), you can get a chart like the one in Fig. 7.12. Fig. 7.12. To make this chart understandable, you should leave only the icons for the empirical data and only the line for the calculated dependence. We have fitted (using the method of least squares, i.e., by minimizing the differences between empirical data and the model) a linear relationship (regression line) that describes the growth of female lake frogs with increasing age. It seems that the obtained distribution has one interesting feature: older frogs are smaller than expected according to the obtained relationship, and middle-aged frogs are either larger or smaller. How to be sure of this? 7.3. Example of analyzing results using conditional formatting. Let's build the distribution of frogs in each age group relative to the size we consider characteristic for that group (according to the model obtained in Fig. 7.11). In column B, we have the actual frog sizes, and in column C, the expected sizes according to our model (~ L=21.9+105.8*A). Let's calculate the relative size of the frog (compared to the expected) in column E. In column F, let's determine which size group each frog belongs to. To do this, we will divide the relative size scale into 6 equal intervals between the minimum and maximum relative sizes. Let's round the relative size of each frog to one of these classes. This will give 7 size classes (less than the first; first; second, third, fourth, fifth, and sixth). The class "less than the first" will be represented by zero; to make all class numbers correspond to natural numbers, we will add one to their numbers. We will get 7 size classes, as shown in Fig. 7.13. By the way, do you understand the reasons why the 1st and 7th classes correspond to twice narrower intervals than the others? These are the same reasons why the first and last classes of the distribution we considered in Fig. 7.1 were twice as small. When constructing distributions for each age (in the lower screenshot in Fig. 7.13, they are in columns I – P), you should specify the row ranges corresponding to individuals of a certain age. Fig. 7.13. The three screenshot fragments show the formulas used to construct the distributions. Figure out how these formulas work! Did you understand why the segment between the minimum and maximum relative size had to be divided into 6 segments (or another even number)? Division into an even number of segments will lead (with standard rounding) to an odd number of classes, one greater. The middle class (in our case, after applying the correction of 1, it is 4) corresponds to individuals whose sizes are approximately as predicted by the model. How to improve the visualization of the distribution in the table we built in columns I – P? Using conditional formatting. To access the option shown in Fig. 7.14, go to "Format / Conditional Formatting... / Color Scale...". Fig. 7.14. Color scale settings dialog. To highlight the 4th size class, which corresponds to medium-sized individuals. Let's put a border around it using "Format / Cells... / Border" or by right-clicking "Format Cells / Border". Ultimately, we get the result shown in Fig. 7.15. Fig. 7.15. Compare this data visualization method with the one shown in Fig. 7.12! These two methods are complementary... In any case, it should be noted that the empirically obtained distribution of individuals by age and size is not very convincingly described by a single linear regression equation... 7.4. Two linear relationships for two intrapopulation strategies Perhaps if one line is not enough, we can try two? This is not the only way to improve the model we used; for example, one could try to complicate the relationship that approximates the existing distribution. However, the author of this course has made such attempts with similar material and was not satisfied with the result. We will consider an approach that corresponds to the concept of intrapopulation ontogenetic strategies, which is presented in these two articles: Shabanov D. A., Korshunov A. V., Kravchenko M. A., Meleshko E. V., Shabanova A. V., Usova E. E. Intrapopulation ontogenetic strategies of precocity and slow growth: definition on the example of anurans // Vestnik of V. N. Karazin Kharkiv National University, series "Biology". — 2014. — Issue 22. — pp. 115-124. Usova E.E., Kravchenko M.A., Shabanov D.A. Intrapopulation ontogenetic strategies in green frogs ( Pelophylax esculentus complex) // Visnik of V. N. Karazin Kharkiv National University, series "Biology". — 2015. — Issue 25. — pp. 223-238. Let's define (next to the previous data block for easier comparison) two model relationships for ~S L and ~B L (from S — small, B — big). Let's calculate the distance from each point to the nearest line (the formula can be seen in Fig. 7.16). We have confirmed that the solution offered by NLPSolver depends on the initial point. Therefore, let's set the initial values for the two lines so that they are close to the previous solution, but one goes below and the other goes above. Fig. 7.16. Model for two lines, ~S L= S a+ S b*A and ~B L= B a+ B b*A. Let's run the "Solver"; set the target and variable cells; select "DEPS Evolutionary Algorithm" or "SCO Evolutionary Algorithm" (Fig. 7.17). You can also try working with the "LibreOffice Swarm Nonlinear Solver (Experimental)" option, but the author of this manual found that the swarm solver produced inadequate solutions when applied to the described problem (it increased, rather than decreased, the target cell value). Fig. 7.17. Conditions are set... Let's run the "Solver"; the algorithm will show how it searches for a solution (Fig. 7.18). Fig. 7.18. "DEPS Evolutionary Algorithm" demonstrates how it searches for a solution. The first search conducted by the author yielded a "suspicious" solution: S b=100. We already know that "DEPS Evolutionary Algorithm" and "SCO Evolutionary Algorithm" change initial values by no more than 4 times; let's run the search again. The final solution is shown in Fig. 7.19, which also includes a graph to visualize this solution. Fig. 7.19. This solution turned out to be optimal after a series of optimization algorithm runs with different initial values. The graph shown in Fig. 7.19 has one significant difference from the similar graph in Fig. 7.12: the lines corresponding to the obtained regression lines are extended to zero age. To do this, three more rows that do not correspond to empirical data have been added to the data range used to build the graph (Fig. 7.10). Fig. 7.20. Rows 109-111 contain points that allow extending the model curves for frog ages not represented in the studied sample. Thus, a model with two regression lines has been built. It undoubtedly describes the diversity of empirical data better than a model with one line. However, perhaps the available data can be described even better. How – the question remains open...