\documentclass[../book.tex]{subfiles}
\graphicspath{{\subfix{.}}}
\begin{document}

\chapter{Fitting Equations to Data}
\label{04EqFit}

\setcounter{section}{-1}  % Start subsections from 0

\input{prereqs/04EqFit_prereqs}

\newpage
\section{Problem}
\setcounter{subsection}{-1}  % Start subsections from 0

The next thing your advisor wants from the UV-vis data are the extinction coefficients for the MLCT bands of the two endpoints (the deprotonated and protonated states), since they will need to be reported in a manuscript. You do a simple experiment in which you measure the absorption spectrum for a series of concentrations, and you know you need to fit this to a line, where the slope allows you to calculate the extinction coefficient. Basically, this comes from Beer's law:
\begin{align*}
    A = \epsilon{c}{l}
\end{align*}
where $A$ is the observed absorption at some specific wavelength, $\epsilon$ is the extinction coefficient, $c$ is the concentration, and $l$ is the path length.

You work with an undergraduate you recently started mentoring, and teach them to perform the data collection on a series of concentrations of your compound.  Later, they send you this data as a series of individual csv files (one for each concentration) and also include an image of the spectrum, taken from the spectrometer. You could do this analysis and fitting in Excel, but if you could do it in Python, it would save time over the Excel workflow that requires opening each file manually, making a chart, adding a trendline, adding the trendline to the display, and using a special function to determine the uncertainty of the fit parameters.  Additionally, there are aspects of fitting in Excel that you don't like. For instance, though Excel will give you the slope and intercept of the best fit line, it doesn't provide estimates of uncertainties in these parameters, and it doesn't show you the residuals of the fit. 


\begin{wrapfigure}{r}{0.5\textwidth}
\vspace{-1em}
    \begin{mdframed}
        \textbf{Why lmfit?}
        
        \noindent A very common practice in experimental science is to try to fit a quantitative model to your data using least squares regression. The Python library SciPy provides functions for this, but they can be a bit tedious to use.  Additionally, analyzing the results (such as uncertainty in estimated parameters) can take some work. The library lmfit streamlines this workflow, and so provides a much nicer interface to work with.  While we think lmfit is the best, it might be worth looking at the library `symfit' as well. 
    \end{mdframed}
\end{wrapfigure}

You talk to \first about this, who says that you absolutely can fit a line to data in Python. They say that there are several ways you can do this, but they think the best is to use a \cite{library} called `lmfit.' ``This library is powerful, and is going to allow you to fit things that are much more complex than a line. It also helps keep track of things like parameter values, uncertainties, residuals, etc. So, it is a good library to know.''


\First says that all you need is the data you already have, and they can even teach you how to quickly select only the csv files and not the images from a folder.  Combined with just a few new ideas about coding, this knowledge can be used to produce a nice result.



\section{Solution}
\setcounter{subsection}{-1}  % Start subsections from 0


\subsection{Planning}

Since there are only two parameters you need to fit (slope and the intercept), you could simply report them (and their uncertainties) as text. But it would be nice to produce a figure that shows your data, the fit to the data, the residuals of the fit. You like the idea of showing the residual as a second plot under the first.

``This is going to require that we understand how \lstinline{lmfit} works and how to include multiple subplots in Plotly. Also, \lstinline{lmfit} has a lot of capabilities, and we are going to be unable to explain everything about it the first time you encounter it, but we will use more of them soon.  Finally, we can introduce a new way top introduce values to code and to identify multiple files to work with."

% Obtaining user input at command line
\subsection{Prompting the user for information using \texttt{input()}}
\label{Prompting_the_user_for_information_using_input}
%\ref{Using input() to get information from the user via the command line}

%\label{Sec:N:04EqFit:input}

``There are going to be many times that you want to run the same basic process, but change only a few specifics within in it.  For instance, we are going to be using absorbance values at the peak of a specific absorption.  But what if you wanted the extinction coefficient for a different spectral feature? You would do the exact same process, but only change the wavelength of interest. Up until now, the only way you could do this is to make a \cite{variable} and assign it the value of interest in the code.  If you wanted to change the value, you would need to change the code.  This is a solution known as \cite{hard_coding}."

``Hard-coding makes it inefficient to change your code dynamically.  A better solution might be to write a program that doesn't know the value of interest ahead of time, but instead asks the user to supply this value only when the code is run. There are several ways you might do this, but one simple way is to use the \lstinline{input()} \cite{function}. The default behavior of \lstinline{input()} is to expect the user to supply an input to the console using their keyboard.  Once the user presses return, the input from the user is \cite{return}ed as a string by the \lstinline{input()} function, which can then be assigned to a variable."

``It is also common to provide \lstinline{input()} a string (as a \cite{positional_argument}) that prompts the user.  For instance, if we wanted to ask what wavelength is going to be of interest for our program we could use:"

\begin{pyin}
l_max = input("What is the position of the feature of interest? ")
\end{pyin}
``If we run this code, we obtain the following output:"

\begin{pyout}
What is the position of the feature of interest? 
\end{pyout}

You notice there is a flashing cursor after the question mark. ``Now, Python is waiting for your input," \first says.  ``If you want to specify a value, you might supply a number like 580.  If you do so, and press return, then that value will be assigned to the variable \lstinline{l_max}."  

You type 580, then hit return.  ``Also, notice that we left a space at the end of the message so that your input would not start immediately after it.  We can verify that it worked by checking the value of \lstinline{l_max} now:"
\begin{pyin}
l_max
\end{pyin}
\begin{pyout}
'580'
\end{pyout}
``One \textit{extremely} important thing to note here is that the value of \lstinline{l_max} is a \textit{string}---as indicated by the fact that the number 580 is placed between single equation marks. You may find this odd, since you input a number.  But remember this is accepting input from a keyboard, and you could just as easily typed in `house,' instead of '580.'  Everything that \lstinline{input()} returns is going to be a string. In our case, it is would be much more convenient if we have a number, and so we are going to need to convert our string '580' into a number. But to do that correctly, we need to understand that, in Python, there are many different kinds of numbers."


% Types of numbers
\subsection{Turning numeral characters into numbers}
\label{Turning_numeral_characters_into_numbers}

``So far, we have treated all numbers as the same. However, Python actually distinguishes between three different kinds of numbers: \citett{int}, \citett{float}, and \citett{complex}. The \lstinline{int} types are integers (numbers without decimal places), the \lstinline{float} types are numbers with decimal places, and \lstinline{complex} types are complex numbers.\footnote{If you want to know more about these numeric data types, you can read about them in Chapter \ref{DataTypes}}  Complex numbers are indicated by using the construction <digits> + <digits>j,\footnote{Many readers may be more familiar with $i = \sqrt{-1}$. Python uses the notation that $j = \sqrt{-1}$ instead.  The reason is that i typically is used for other purposes in programming, such as keeping track of index values.} where the first set of digits are the real portion, and the second are the complex portion. We have already been using \lstinline{int}s and \lstinline{float}s this whole time, though we have not had to worry about it until now."

``For now, we simply need to know two things.  First, we can create an \cite{instance} of each type using the function with the same name: \lstinline{int()} will create integers, \lstinline{float()} will create floats, and \lstinline{complex()} will create complex numbers.  You've seen that normally, we don't have to specify the type, and Python just figures it out. However, when converting between data types (say, from a \cite{string} to a number type), then we have to worry about it. Second, we need to know that, from the perspective of Python, these numbers are different.  When we want to do something like convert a string to a number, we need to specify which type of number we want.  For instance, if we had the string "2.6" and wanted to convert it, we will get different results for each number type:"
\begin{pyin}
textnumber = "2.6" # this is a string
\end{pyin}
\begin{pyin}
int(textnumber)
\end{pyin}
\begin{pyout}
2
\end{pyout}
\begin{pyin}
float(textnumber)
\end{pyin}
\begin{pyout}
2.6
\end{pyout}
\begin{pyin}
complex(textnumber)
\end{pyin}
\begin{pyout}
(2.6+0j)
\end{pyout}
``If we really want to have numbers with decimals, we need to ensure we are using \lstinline{float()} to create them."

\subsection{Finding the element closest to a given number in a Numpy array}
\label{Finding_the_element_closest_to_a_given_number_in_a_Numpy_array}

The next thing you start to worry about is, what happens if the user specifies a wavelength that isn't exactly one of the data points in the spectra. You mention this to \first, who says, ``There are going to be many times that we want to find the index of the value that is the closest to a target value.  A simple way to do this for \cite{Numpy array}s is to subtract the target value from the whole array, take the absolute value of the differences, and then find the minimum of that.  The number that was closest to the target value will result in the smallest number." 

``Because we know we can use this approach, the Numpy \cite{function} \lstinline{argmin()} becomes very useful.  What this does is examine a Numpy array and return the \cite{index} for the \textit{first} number that is the smallest in the array.  We can use this index to then get the number of interest."

\begin{lstlisting}
import numpy as np

arr = np.array([0, 3, 1, 2.3, 3.1])
target = 2.7
ind = np.argmin(abs(arr - target))
print("closest to", target, "is", arr[ind], "at index", ind)
\end{lstlisting}
\begin{pyout}
closest to 2.7 is 3.0 at index 1
\end{pyout}
``Here, we get the index for the number closest to 2.7, and then get that number using the index. One important thing to realize about \lstinline{argmin()} is that it does return only the first instance of this smallest value. So, if we do:"
\begin{lstlisting}
arr = np.array([0, 3, 1, 3, 3.1])
target = 2.7
ind = np.argmin(abs(arr - target))
print("closest to", target, "is", arr[ind], "at index", ind)
\end{lstlisting}
``we will obtain the exact same result as before, even though these is also a value of `3' at index 3."
\begin{pyout}
closest to 2.7 is 3.0 at index 1
\end{pyout}
``This is worth keeping in mind as you write your codes."

``If there is an \lstinline{np.argmin()}, it should not be surprising that there is a \lstinline{np.argmax()} as well. This does the same thing, but returns the first occurance of the index for the maximum value in the array."


%f-strings
\subsection{Inserting calculation results into strings using f-strings}
\label{Inserting_calculation_results_into_strings_using_f-strings}

``It is going to be convenient to write out the results of our fit using strings. But since we also have uncertainty in our fit we would like to include that. Ideally, would be nice to write a sentence like "<value> $\pm$ <uncertainty>". In our program, <center> and <uncertainty> will be held in \cite{variable}s. We do have a way to print this, using a construction like \lstinline{print(<center value>, " +/- ", <uncertainty>)}, but there are two problems with this solution.  First, it is a bit cumbersome to write and read.  Second, this doesn't actually produce a \cite{string}, so this construction provides no way to pass the sentence we want to any other function, like Plotly's \lstinline{add_annotation()}."

``There are several ways to solve the problem of getting a string object with values from variables,'' \First says, ``But I think the best option is using \cite{f_string}s.  These are strings, but we place an \lstinline{f} prior to the first set of quotation marks. We can use this to produce a 'normal' string:"
\begin{pyin}
 f"Hi"
\end{pyin}
\begin{pyout}
'Hi'
\end{pyout}
``But we can also place variable names between curly brackets (\{\}) and the variables will be evaluated." 
\begin{pyin}
pi = 3.14159
\end{pyin}
\begin{pyin}
f"A nice approximation of pi is {pi}."
\end{pyin}
\begin{pyout}
'A nice approximation of pi is 3.14159.'
\end{pyout}
``Actually, f-strings do more than this, since really what is happening is that any Python \cite{expression} we place between the curly brackets is evaluated and the value included in the string.''

\begin{pyin}
f'The square of pi is {pi**2}'
\end{pyin}
\begin{pyout}
'The square of pi is 9.869587728099999'
\end{pyout}
``These evaluations could include function calls, or really any valid Python statement. As you can imagine, f-strings can be useful in a huge range of contexts, even when you're not outputting text." 

% formatting numbers
\subsection{Formatting the representation of numbers in strings}
\label{Formatting_the_representation_of_numbers_in_strings}

``If we are going to be placing numeric values inside of \cite{string}s, then it is going to quickly become important to understand how to format these values.  For instance, imagine you had 10 g of a chemical you wanted to equally divide up between 7 vials. You could use f-strings to provide a nice statement."

\begin{pyin}
f"Place {10/7} grams in each vial."
\end{pyin}
\begin{pyout}
'Place 1.4285714285714286 grams in each vial.'
\end{pyout}
``This is not a great output, because it is hard to read and you likely do not have a balance that reads with attogram precision. Even if you did, there is no need for such precision. Instead, weighing to 1 mg might be sufficient, so we only need to report this to 3 decimal places. We can format the string, by placing a colon (:) after the value and then formatting instructions before the end of the curly brackets."

\begin{pyin}
f"Place {10/7:.3f} grams in each vial."
\end{pyin}
\begin{pyout}
'Place 1.429 grams in each vial.'
\end{pyout}
You note that the number has even been rounded, so that you don't have to do it yourself as a separate step. 

``In this usage, the 0.3 means that you want 3 decimal places, and the 'f' indicates that this should be formatted as a standard \ul{f}loat. There are other options where we can ask for scientific notation, \cite{e_notation}, etc."\footnote{More information can be found in technical Chapter \ref{Text}.}


\subsection{Generating symbols that aren't on your keyboard using \texttt{codechembook.symbols}}
\label{Generating_symbols_that_aren't_on_your_keyboard_using_codechembook.symbols}

``As we do more sophisticated things in Python, we are going to want to produce sophisticated \textit{looking} text. For instance, rather than writing +/- to indicate plus-minus values, we would want to use the character: $\pm$."

You note that your keyboard does not have a $\pm$ key, so you wonder how you could type this character. \First says that what you need is the designation for this character that your computer uses. Since a computer views everything as a number, characters are also just numbers. The exact number depends on the type of \cite{binary_encoding} being used, but a common one (and one that Plotly can use) is called \cite{unicode}.\footnote{If you want a more complete discussion of this, you can refer to Chapter \ref{Text}.} In unicode, all characters are represented by a string of \cite{hexadecimal} numbers that follow \lstinline{U+}. For instance, the character `t' is \lstinline{U+0074}. Thus, when using unicode encoding, the keys on your keyboards are just shortcuts to these designations.  If you want to render a character that is \textit{not} on your keyboard, you could look up what the corresponding unicode is, and then enter it directly. For instance, the unicode for $\pm$ is \lstinline{U+00B1}.

``This approach will work, but can be cumbersome, so the \cite{codechembook} \cite{library} has a module with shortcuts that use variable names as aliases for the unicdoe strings of common glyphs that you might need.  For instance, the plus-minus glyph can be accessed using:"
\begin{pyin}
codechembook.symbols.math.plusminus
\end{pyin}
``and this will output $\pm$."

``There are useful math glyphs in the \lstinline{codechembook.symbols.math}, such as \lstinline{proportional_to} ($\propto$) and \lstinline{equivalent_to} ($\equiv$)."

\begin{mdframed}
    \textbf{Finding and Using Unicode Characters}
    
    \noindent While \lstinline{codechembook.symbols} provides a number of useful symbols, it only a small fraction of all possible symbols (including an ever-growing number of emojis).  To find other symbols, you can use lookup websites like \href{https://unicodelookup.com/}{https://unicodelookup.com/}, which will give you the \cite{hexadecimal} code that specifies the symbol.  Unicode characters can be used in strings by including the hexadecimal code (minus the beginning \lstinline{0x} after \lstinline{\u}: the superscript plus, for instance, can be inserted into a string by placing \lstinline{\u207A} where the superscript plus should be.
\end{mdframed}

``The \lstinline{codechembook.symbols} module has many different symbols, which are organized into several sub-modules.\footnote{You can find a description on the \website and they are also discoverable when using them in an IDE, meaning that auto-complete hints should pop up with them.}  These are \lstinline{math}, \lstinline{greek}, \lstinline{script}, \lstinline{chem}, and \lstinline{typography}. We have already seen examples from \lstinline{math}.  The \lstinline{greek} sub-module contains Greek letters, such as \lstinline{Gamma} ($\Gamma$) and \lstinline{gamma} ($\gamma$). The \lstinline{script} sub-module contains script versions of Latin letters, such as \lstinline{l} ($\ell$) and \lstinline{H} ($\mathcal{H}$). The \lstinline{chem} sub-module contains symbols that might be of use to chemistry, such as \lstinline{degree} (°), \lstinline{angstrom} ($\AA$), and \lstinline{wavenumber} (cm$^{-1}$). Finally, the typography sub-module contains glyphs that do not fall into the above categories, such as \lstinline{copyright_} (©) and \lstinline{ellipsis} (…)."


%formulas
\subsection{Producing formatted chemical formulas using \texttt{codechembook.quickHTMLFormula()}}
\label{Producing_formatted_chemical_formulas_using_texttt_codechembook.quickHTMLFormula()}

``Since we are chemists, we often want to include chemical formulas on plots that we make---for instance, we can add a title to a Plotly plot using the \lstinline{title} \cite{keyword_argument} inside the \lstinline{layout()} \cite{method} of the figure object. We know that, for annotations and labels in Plotly, we need an HTML formatted string.   We already know how to produce a properly formatted chemical formula using the HTML superscript and subscript tags \texttt{<sup></sup>} and \texttt{<sub></sub>} or the quicker functionality of \lstinline{codechembook.symbols.typsettingHTML}. However, doing it this way is tedious because you have to format each \cite{string} element separately. Fortunately, there are tools out there that can parse a chemical formula and produce a useful string for you to include in a plot or anything else. One approach is included in the \lstinline{codechembook.quickTools} module.''

``We can use \lstinline{codechembook.quickTools.quickHMTLFormula()} to produce these HTML-formatted chemical formulas.  To use this function, we simply pass the formula you want as a string for the first \cite{positional_argument}. For instance, if we wanted the condensed formula for $n$-octane, we could type:"

\begin{pyin}
from codechembook.quickTools import quickHTMLFormula as qf
\end{pyin}
\begin{pyin}
qf("CH3(CH)6CH3")
\end{pyin}
\begin{pyout}
'CH<sub>3</sub>(CH)<sub>6</sub>CH<sub>3</sub>'
\end{pyout}
``You can see that this function returns a string that automatically contains HTML code specifying the subscripts needed.  If you want to specify charge, you don't do it in the chemical formula, but use the \lstinline{charge} keyword argument.
\begin{pyin}
qf("SO4", charge = -2)
\end{pyin}
\begin{pyout}
'SO<sub>4</sub><sup>2-</sup>'
\end{pyout}

``If you want strings formatted for other applications, there are other functions. For instance, you can get a unicode-formatted string using \lstinline{quickUnicodeFormula()} if you want to output, say to the command line instead of Plotly.  You can get a \Latex-formatted string using \lstinline{quickLatexFormula()}.  For instance:''
\begin{pyin}
from codechembook.quickTools import quickUnicodeFormula as quf
\end{pyin}
\begin{pyin}
quf("SO4", charge = -2)
\end{pyin}
which will produce: 'SO$_4$²$^-$'

These functions from \lstinline{codechembook} are implementations of the \lstinline{chemformula} \cite{library},\footnote{https://github.com/molshape/ChemFormula} but we have included them in our library, because \lstinline{chemformula} was not available on Conda (at the time of this writing). The \lstinline{ChemFormula} object that the \lstinline{chemformula} library produces also contains a number of other useful attributes and methods.  You can also import the \lstinline{codechembook.chemformula} module, which will give you this functionality, and read the documentation for the \lstinline{chemformula} library on how to use it. 



%glob
\subsection{Finding paths to files in folders using \texttt{glob()}}
\label{Finding_paths_to_files_in_folders_using_glob}

``Sometimes, you will want your code to select files for you, rather than specifying them yourself using code or the file dialog.  This is especially true if you have a folder with many different types of files in it.  For this, the \lstinline{glob()} \cite{method} of the \acite{Path} \cite{object} is very useful.  It accepts an argument that specifies a pattern, and then \lstinline{glob()} looks through the  Path object and finds all items that match this pattern. 

However, we would rather not have to know or specify the entire path for every file, since your files will have unique names.  Instead, we often specify a pattern to match file names to that uses what is called a '\cite{wildcard}.'  Essentially, any file name that contains the text specified in the pattern, \textit{and} that contains text that can be replaced by the wildcard ('*') and still match the pattern, will be valid. For instance, imagine we wanted to find all of the Python script files (files with \cite{extension} '.py') in a folder, named <folder>.  If we make a Path object to this folder, then we can use: \lstinline{<folder>.glob('*.py')}."

``In this usage, \lstinline{glob()} will look through the <folder> object and find all files that have '<anything>.py'.  Since the '.py' is specified, the files must contain it, and since there is a wildcard ('*') before that, we will find all files that end in '.py.'  Since Python files end in '.py,' \lstinline{glob()} will find all the Python files in <folder>."

``One technical aspect to be aware of is that \lstinline{glob()} produces an object called a \cite{generator}. We encountered this for the first time when we introduced the \lstinline{range()} function.\footnote{See Section \ref{Generating_sequences_of_numbers_using_range}.} A generator is an object that doesn't have the objects you will use, until you ask it for them (i.e., it generates them on demand).  One way to get these objects all at once is to ask for them as a \cite{list} using: \lstinline{list(<glob generator object>)}. Here, the generator will generate \cite{string}s that contain all the filenames and their paths, and put them in a list."

``Finally, we can use the \lstinline{codechembook.quickTools} module function \lstinline{quickSelectFolder()} to open a dialog box that allows us to get the path to a folder that we want to glob, or any other scenario where we might want to have the path to a folder instead of file(s)."

% splitting strings
\subsection{Dividing a string into multiple substrings with \texttt{split()}}
\label{Dividing_a_string_into_multiple_substrings_with_split}

Now you start to think about getting the values for the concentration of each spectrum.  You remember that they are recorded in the file names of your spectra.  You could make a new file with these separately, and read that, or you could hard code them.  ``Or!" \first says, ``this is an example of a time where it is useful to break apart text \cite{string}s to get their contents.  As a simple example, think about what we would do with a csv file, where values are separated by commas.  We would break apart the text where these commas are and put the pieces into a list. We can accomplish this using the \lstinline{split()} \cite{method} of a string \cite{object}. This \cite{function} accepts a \cite{positional_argument}, which is the pattern we will split on."

\begin{pyin}
data = "1325, 32512, 232314, 12311,"
\end{pyin}
\begin{pyin}
split_data = data.split(",")
\end{pyin}
\begin{pyin}
print(split_data)
\end{pyin}
\begin{pyout}
["1325", " 32512", " 232314", " 12311", ""]
\end{pyout}
Here you can see that \lstinline{split()} \cite{return}ed a \cite{list} of strings, where the strings are the bits that were separated by commas in the original string. You also notice that the list items are still strings, so they are not numbers (from a \cite{data_type} perspective).  Finally, it is clear that the text split \textit{exactly} on the supplied argument.  Since you provided "," with no spaces, it has split the string only on the comma, and preserved the spaces after the commas. At first, it is confusing that list ends with an empty string.  You realize that the original string ended with a comma, so there was nothing that followed it, and the split function returns the "nothing" that followed it, but as a string---in other words the empty string. 

``In your case, all you would do is split your file name using whatever character is in front of the concentration!"





% \textcolor{red}{IT APPEARS THAT WE DON'T INTRODUCE CONDITIONAL EXPRESSIONS TIL NEXT CHAPTER}

% There are also times that we wish to find values that are neither the minimum nor the maximum. There are several routines one might use, but we will introduce only one: \lstinline{np.where()}.

% As we wish to use it, \lstinline{np.where()} accepts one positional argument: a logical condition involving an array.  So, imagine we had an array that included values from 0 to 10 and wanted to find the index of all values that were greater than 2:

% \begin{pyin}
% array2 = np.array(range(0,10))
% \end{pyin}
% \begin{pyin}
% np.where(array2 > 2)
% \end{pyin}
% \begin{pyout}
% (array([3, 4, 5, 6, 7, 8, 9]),)
% \end{pyout}
% You can see this returns a tuple that has two entries: an array with the indices of the elements that met the condition, and an empty element.  Thus, if you only wanted the indices, you would use:

% \begin{pyin}
% np.where(array2 > 2)[0]
% \end{pyin}
% To get the first element of the tuple, which is the array of indices. 







%lmfit
\subsection{Fitting models to data with \texttt{lmfit}}
\label{Fitting_models_to_data_with_lmfit}

You recall that the process of fitting \cite{model}s to data includes the following steps: identify a model that produces our expected behavior, determine the fitting parameters that control the behavior of our model, and optimize the values for these parameters until the predicted behavior is `close enough'\footnote{There are multiple ways that this is defined. The default for lmfit is that if the difference in the sum of residuals squared between two parameter guesses is below some threshold. If you wish to have more information, consult the online documentation for lmfit.} to the experimental behavior.  At this point the fitting is finished and you use the final parameter values as the best-fitting values for those parameters.  In an ideal world, you can also analyze the process of the fitting in order to estimate uncertainties in these final parameter values, though you remember this being mathematically complex. 

\First agrees.  ``Fitting can be quite involved, but fortunately, there is a \cite{library} that can simplify this for us: \citett{lmfit}. If you really want to understand how lmfit works, then you can read the online documentation.  However, the basics of its use are as follows:
\begin{enumerate}
    \item Create the model we wish to use, and the parameters needed for that model.  In our case, we want to use a linear response model, and so the parameters are the slope and intercept of a line. 
    
    \item Define the parameters and make an initial guess at their values. 

    \item Fit the model object to our experimental data, using these initial parameter value guesses. This will generate a `ModelResult' \cite{object} that contains a significant amount of information about the fit. We can use \cite{attribute}s and \cite{method}s of this object to access this information.
\end{enumerate}
Again, more information can be found in the online documentation for lmfit. However, this list provides enough structure to follow how we will use it below." 



\subsection{Analyzing fit results in \texttt{lmfit}}
\label{Analyzing_fit_results_in_lmfit}

``One nice thing about \lstinline{lmfit} is that the result of the fit---a ModelResult \cite{object}---has a \cite{method} for generating a report of the fit: \lstinline{fit_report()}. Using this method without any \cite{argument}s will produce a formatted report that gives important information regarding the fit.  Specifically, it will provide: the model that was fit to the data, the method of fitting, several goodness of fit parameters (including reduced-$\chi^2$ and R$^2$), the final value of the parameters and their uncertainties, and any correlation between the parameters. Thus, simply by invoking this function, we can see most of the information we would be concerned about."

``Though metrics like reduced-$\chi^2$ and R$^2$ are commonly used to indicate goodness of fit, there are other considerations that are just as important.  For instance: does the resulting fit line make sense? Is there a pattern in the residuals (the difference between the fit line and the data) that needs to be explained? For this reason, we recommend \textit{at a minimum}, always plotting the data, the fit result, and the residual.  Fortunately, \lstinline{lmfit} also provides ready means to do so. The fit result object has many \cite{attribute}s, but the ones we will care most about are \lstinline{best_fit} and \lstinline{residual}. They provide the $y$-values for the best fit and the residuals of the fit, respectively."


% subplots
\subsection{Building figures with multiple subplots in Plotly}
\label{Building_figures_with_multiple_subplots_in_Plotly}

``We are already using the \lstinline{make_subplots()} function. So far, we have not provided it any arguments, so it produces the default of a single subplot.  However, we can tell it how many subplots we want to use, which is done using the \lstinline{rows} and \lstinline{cols} \cite{keyword_argument}s (notice the `\lstinline{s}' at the end of these keywords). Each of these keywords accepts an integer that specifies how many rows and columns of subplots you want. Then, when using other functions that modify a plot (e.g., \lstinline{add_trace()}, \lstinline{update_xaxes()}, etc.) we use the \lstinline{row} and \lstinline{col} keyword arguments (notice the \textit{lack} of an `\lstinline{s}' at the end of these keywords) to specify which subplot to apply them to.  Counting for the subplots starts at 1, row counting starts at the top, and column counting starts on the left. Thus, \lstinline{row=1} and \lstinline{col=1} always selects the top left subplot."\footnote{Much more information on this can be found in Chapter \ref{ch:plotly}.}


% \subsection{It is good to understand the structure of your data}

% "Alright, the last thing we need to do is understand the structure of our data, because that will determine how we need to process the files to get the data out."

% \First opens up one of the files, and it looks like the following:
% \begin{figure}[h]
%     \centering
%     \includegraphics[width = 1.3in]{ChapNarrative/04EqFit/figures/EqFitData.png}
%     \caption{Screenshot of part of the csv data that will be plotted, opened up in a text editor. The relevant points here are that there is a \cite{header} in the first row of the file, and then columns containing \cite{integer} values for wavelength, and absorbance values using numbers with decimals. }
%     \label{fig:EqFit:csv}
% \end{figure}

% \\"So you can see that you have a csv file with a \cite{header} in the first row, the wavelengths are in the first column, followed directly by a comma, then there are the absorbances, with no comma following them. We can work with this structure to produce a code that reads the data, fits a linear model to the data, prints the report of this fit, and the plots the fit and its residuals. So let's do it!"



% CODE
\section{Code}
You are now ready to produce code that can ask the user for the wavelength they wish to find the extinction coefficient for, read the data from csv files and extract the targeted data, and then fit this to a linear model, outputting the results, and a plot of the results.

\setcounter{subsection}{-1}  % Start subsections from 0

\begin{figure}[H]
    \centering
    \includegraphics[]{ChapNarrative/04EqFit/figures/FlowChartCh5.pdf}
    \caption{Flowchart for the code to determine the extinction coefficient of your catalyst.  It reads the csv files in the Extinction data folder containing UV/vis spectra at a range of concentrations, extracts the absorbance value at a specified wavelength, and fits this data to a linear model. }
    \label{fig:EqFit:ActivityDiagram}
\end{figure}

\lstinputlisting[language=Python]{../../Distribution/Code by Chapter/Ch05 - Fitting Equations to Data/ExtinctionCalc.py}

You run this code, navigate to the folder ``Titration," which is a subfolder of ``Catalyst," and select it.  You enter 581 when the code asks you for the wavelength at which to calibrate, and it quickly produces the following fit report output and figure:

\begin{figure}[h]
    \centering
    \raisebox{4em}{\includegraphics[width = 2.4in]{ChapNarrative/04EqFit/figures/eqFitReport.png}}
  %\hspace*{.2in}
  \raisebox{-0em}{\includegraphics[width = 3.4in]{ChapNarrative/04EqFit/figures/eqFitPlot.png}}    \caption{Results obtained using the lmfit \cite{library} to fit a linear model to data of absorbance versus concentration. The left shows the output from the \lstinline{.fit_report()} method of the lmfit result \cite{object}, as output into a \cite{console}. The right shows a plot of the data, fit, and residuals of the fit, plotted in two subplots using Plotly.}
    \label{fig:EqFit:FitResult}
\end{figure}


%\raisebox{-0.5\height}


\subsection{Line-by-line explanation}

\begin{footnotesize}
\begin{itemize}
    \item [9--13:] Import the libraries we need.  One new library is \lstinline{lmfit}, from which we import the \lstinline{LinearModel} from the \lstinline{models} \ul{class}. This will provide the behavior we want to fit to our data.  We also introduce the \lstinline{quickSelectFolder()} function from the 
    \lstinline{codechembook} library. This function allows the user to select a folder using the native file explorer, and returns the path to this folder as a Pathlib Path object. We further import the \lstinline{quickHTMLFormula()} function to handle formatting of our chemical formula for the title of the plot.
    Finally, we import the \lstinline{symbols} module of \lstinline{codechembook}, which we are going to use to get \cite{unicode} \cite{string} for characters that are commonly used in science.  

    \item [16:] We are going to have the user identify a folder.  Since we have been using the console to communicate with the user, we can print a message here that tells them what we want them to be doing.  This can be helpful when people other than the coder run the code, or if you come back to your code after a long time. Moving forward, we will explore other means to prompt the user, but for now this is a nice solution.

    \item [17:] We use the \lstinline{quickSelectFolder()} function to have the user identify a folder where the data is.  This will return a Path object, so we can use the \lstinline{glob()} method of this Path object. This method will find all files that match the pattern supplied.  In this case, the pattern is a wildcard (\lstinline{*}) followed by \lstinline{.csv}, so \lstinline{glob()} will match to all files ending in ".csv". Glob returns a \cite{generator} for the specified Path objects in this folder, and we pass this result to \lstinline{sorted} function, so that this list is sorted, and then we assign the result to the variable \lstinline{filenames}. It is not strictly necessary that the files be in order, but it does make things a bit nicer down the line and costs very little time or typing. 

    \item [20:] We are \cite{nesting} functions here. By this we mean that we are taking one function (\lstinline{input()}) and then using the value that is returned as the argument for another function (\lstinline{float()}). Nested functions evaluate from the innermost set of parentheses out.  In this case, the effect is to obtain input from a user, which we expect to be a number, and then convert this input from a string (the \cite{data_type} \lstinline{input()} returns) to a \cite{float}. If the user supplies a data type that cannot be converted to a float, then the code will stop and automatically generate an \cite{error} message.

    When running this code on your own, we think a good wavelength to use here is 581 nm, but feel free to play around with other values!

    \item [23:] To form a calibration curve, we need two sets of data: the concentrations of our species and the absorption at these concentrations. We make two empty lists to which we will add these values, as we pull them from files and assign them to two variables on the same line, using \cite{unpack}ing. 

    \item [24--29:] We start a \citett{for} \cite{loop}, in which we go through each item in our list of filenames one by one.  This means we are iterating through the Path objects that are the paths to our data files. 
    \begin{itemize}

        \item[25:] We need to know the concentration of the species.  This is contained in our filenames, so we can extract it from the filenames. Because the \lstinline{f} variable points to a Path object, we can get just the name (without the extension) using the \lstinline{stem} \ul{attribute} of the Path object. Once we have the filename, we can use the \lstinline{split()} method of the string to split it into a list where each item in the list is a part of the filename before or after a ``\_".  We then get the \textit{second} item in the list (index 1) returned from the \lstinline{split()} method, which is the concentration.  However, this is returned as a string, and we want it to be a float, so all the prior commands (which return the concentration as a string) are contained in the \lstinline{float()} function to convert the string to a floating point number. Finally, all of this is contained in the \lstinline{append()} method of our \lstinline{conc} list.  The end result of this line is that we take the file we are currently on in the \lstinline{for} loop, and we extract the concentration as a float and add it to our list of concentration values. It is worth noting that this only works because of the file naming convention that we have used---so in your own work, it might be worth carefully considering how you name your files. 

        \item[26:] We use the \lstinline{genfromtxt()} function, just as we have in past chapters, to extract the data of interest from the file, unpacking the returned arrays into two variables. These hold the $x$ and $y$ data.

        \item [27:] We use the \lstinline{argmin()} function of Numpy to find the index corresponding to the $x$-value that is closest to the value of \lstinline{l_max}. The way this works is that we subtract the value of \lstinline{l_max} from our Numpy array of $x$-values.  We then take the absolute value of this result, using the \lstinline{abs()} function.  In the resulting Numpy array, the item with the smallest value will be the item at the same index as the item in the original array that was closest in value to \lstinline{l_max}.  We then use \lstinline{argmin()} to find the index of this smallest value, and assign this index to a new variable. 

        \item [28:] Using the index we obtained in the line above, we obtain the value held in the Numpy array of $y$-values at that same index. This results in pulling the $y$-value for the index that was closest to our target \lstinline{l_max}.  We add this to the \lstinline{absorb} list using the \lstinline{append()} \ul{method} of the list. 
    \end{itemize}
    
    \item [31:] We create a lmfit Model object, using the \lstinline{LinearModel()} function that we imported on line 10. This Model includes instructions to create a line and manage the parameters \textit{slope} and \textit{intercept}. We assign this to a variable for use later on. 

    \item [32:] We invoke the \lstinline{guess()} method of the \lstinline{LinearModel} object, passing it the $y$-values as a positional argument and $x$-values as a keyword argument.  This method then determines initial guesses for the model parameters.  Since this is a linear model, the parameters are the slope and intercept of the line. We assign these initial guesses object to a variable. 

    \item [33:] We use the \lstinline{fit()} method of the \lstinline{LinearModel} object to fit the model to our data.  So we provide the $y$-values and the initial values of the parameters as positional arguments, while we also provide the $x$-values as a keyword argument. This will return a fit object, which we assign to a variable.

    \item [34:] Using the \lstinline{fit_report()} method of the result object, we get the fit report, and then print it to the console. This produces a report from which we can read the goodness of fit (R$^2>0.99$ in this case), as well as the final values and uncertainties for the parameters. 

    \item [37:] Print out the molar absorptivity to the user in the console.  Here we make use of f-strings and the \lstinline{chemSymbols} feature of the \lstinline{codechembook} library to give us the $\pm$ symbol and the superscripts we will use on the units.

    \item [40:] We make a Plotly Figure object that will hold 2 subplots. We want 2 rows and 1 column, which we specify using keyword arguments. This will produce a vertical arrangement of the two subplots. 

    \item [43--45:] We add traces to the Figure object. The first two are added to the top plot (row = 1, col = 1), while the third is added to the lower plot (row = 2, col = 1). All three use the same $x$-values, which we can access using the unintuitive \lstinline{result.userkws} attribute.  This attribute is a dictionary in which the independent variable array is stored using the key that is specified for the fit model, in this case, simply \lstinline{'x'}.  The first line adds the results of our fit (using the \lstinline{best_fit} \cite{attribute} of the result object), the second adds the experimental $y$-values, and the last adds the residuals (using the \lstinline{residuals} attribute of the result object. We also choose to show the experimental data using markers only, the fit result using lines only, and the residuals using markers and lines. Since the last of these is the default behavior, we do not need to specify a mode.  Finally, we include none of these traces in the legend. 

    \item [48--51:] We are going to annotate our plot with three pieces of information: the final slope value (and uncertainty), the final intercept value (and uncertainty), and the r-squared value of the fit. There are many ways to do this, but the simplest is going to be to use a multi-line f-string. Multi-line strings are placed between """ characters, and this allows us to specify strings that span multiple lines, which will help us visualize the string we are writing.  Normally, the line breaks in the multi-line strings also translate to breaks in the string when printed or otherwise represented. However, because we are going to put this string as an annotation in a Plotly figure, `normal' line breaks do not produce line breaks.  Instead, we have to the use HMTL tags for breaks: <br>. 

    We will create a line for three pieces of information from the fit: the slope (and uncertainty), the intercept (and uncertainty), and the r-squared value. These are passed to the f-string, put between curly brackets (\{\}). We can use the \lstinline{.value} attribute of the parameters accessed via a dictionary in the fit result object to access the values, and the \lstinline{.stderr} attribute to access the uncertainty in these parameter values. 

    When we supply these values, we can also format them, providing the information after a `:' within the curly brackets, after the value we are formatting. In this case, we specify we want \cite{e_notation} with two decimals of precision using `.2e'.

    One final note: we also use the chemical symbols function from the \lstinline{codechembook} library to provide quick access to the $\pm$ symbol in our plot. 

    We assign this string to a variable, so that we can use it below. 

    \item [53--56:] We add the above string to the Plotly figure using the \lstinline{add_annotation()} method of the Figure object. We provide the string we just created above, and specify that we do not want an arrow and where to place this annotation on the plot, using the $x$ and $y$ axis coordinates.  

    \item [59--62:] Also using the \lstinline{add_annotation()} method, we add the molar absorptivity to the plot. This time, we construct the string we want inside the function call.

    \item [65--67] We add a title to the $x$ and $y$ axes. However, the title for the $y$-axis is different for the top (experimental data) and bottom (residuals), so we specify a title for each $y$ axis individually, using the \lstinline{row} and \lstinline{col} keyword arguments.  Note that, when assigning the title to the $x$-axis, we do not specify a row or column, and so the title is applied to \textit{all} of the $x$-axes in the plot. 

    \item [68:] We apply the \lstinline{simple_white} template to our plot. We also add a title using the \lstinline{title} keyword argument and supply a properly formatted HTML string, returned by the \lstinline{quickHTMLFormula()} function.

    \item [69:] We show the result of the plot. 
\end{itemize}
\end{footnotesize}


\hrule
\input{Exercises/Chapter 05 Fitting equations to data}
\hrule

\end{document}
