path: root/matlab_usage_source/index.t2t

Matlab usage - short (may be too short) introduction
by EugeniyMikhailov



%!includeconf: ../configs/config.t2t
%!style(html): ../evmik.css

= Why Matlab =

At the very least, it is a great replacement of a  hand held calculator. Once
you need  to process multiple  data point it  became tedious with  a simple
calculator. Once  you master Matlab  it will not matter  if it is  one data
point or millions of them. It  can do complex math, plotting, fitting, etc.
It is also the complete programming language which can stand by itself.
Also, Matlab has **excellent help documentation** and a lot of tutorials at the web.

There are free alternatives. [Octave  http://www.gnu.org/software/octave/]
is  one example. This what I personally use
instead of Matlab all the time. But it is a bit difficult with installation
and it  might be  scary to  see no  GUI and  only a  command prompt  for a
beginner.

= Getting Matlab =

Please visit W&M IT [software support page http://www.wm.edu/offices/it/a-z/software/]
and download Matlab from appropriate "Licensed Software >> Math & Statistics
Software" section. They have several available versions. Either one is fine.
Since we are learning Matlab, we will not  have time to go to fancy toolboxes
which Matlab provides/removes with new releases. 


= Statistics functions =

== Mean of a set ==
Often you will need to find the mean and standard deviation of a particular
sample set. Let's say we estimated  a resistor value (**R**) with different
method and obtained the following experimental values.

| **R**  | 1.20 | 1.10 | 1.11 | 1.05 | 1.05 | 1.15 |

First we transfer it to a Matlab vector variable
``` R=[ 1.20, 1.10, 1.11, 1.05, 1.05, 1.15 ]

Notice that for a vector variable I use **,** as a separator.

Now we ready to calculate mean of this set
``` mean(R)
with answer 1.1100

== Standard deviation of a set ==
Standard deviation normalized with **N-1** samples  (so called unbiased estimator)
is done with
``` std(R)
with answer 0.058310. If you need normalization with **N** samples do
``` std(R,1)
with answer 0.053229

In both  cases there is  to many  significant figures which  is unjustified
from physics point of view. So remember, Matlab just does the math for you,
the analysis and presentation of results is on you.



= How to make a simple plot =
Suppose during the measurements we obtain the following data
for each current value **I** we measured voltage **V**
with some uncertainty in voltage measurement **dV**.
Our experimental data is represented in the following table.

|| I   | V   | dV  |
|  0.1 | 1.0 | 0.2 |
|  1.0 | 2.2 | 0.1 |
|  1.8 | 3.2 | 0.2 |
|  3.3 | 3.8 | 0.2 |
|  4.0 | 4.9 | 0.3 |

First of all,  we need to transfer  it to the Matlab column vector variables. Pay
attention to the delimiter **;** between numbers

``` 
I  = [  .1; 1.0; 1.8; 3.3; 4.0 ]
V  = [ 1.0; 2.2; 3.2; 3.8; 4.9 ]
dV = [ 0.2; 0.1; 0.2; 0.2; 0.3 ]
```

Now we make a simple plot without error bars to get  the basic, text after
**%** is considered as a comment and Matlab will disregard it.

```
figure(1) % all plotting output will go to the particular window
plot(I,V, 'x') 
```

[simple_unlabeled_plot.png]

Easy,  is not  it? Here  were mark  our data  point with  crosses and  used
**'''x'''** parameter to specify this. You can do other markers as well.


But  **every plot  must be  present with  title, proper  axes labels  and units**.
Note **'** for the string specifiers

``` 
title('Dependence of voltage on current')
xlabel('Current (A)')
ylabel('Voltage (V)')
``` 

[simple_labeled_plot.png]

No our plot is nice and pretty.


= Fixing the range of the axes =

Have a look at above plot. We see that one of the axes (y) does not contain
origin point  (y=0). Such plots are  very common in finances  and design to
fool a reader, but they are  generally not acceptable in scientific reports
and  definitely not  in  my class  unless  you  have a  good  reason to  do
otherwise.

The fix is very easy. Note that we are actually sending a vector [0,5] to the
**``ylim``** function (fir x-axis we need to use **``xlim``** function).
``` ylim([0,5])

Now our plot is good to present

[simple_labeled_and_proper_limits_plot.png]


= Plotting with error bars =
Let's open a new window for the plot

``` figure(2)

Now we make plot with error bars assuming that they are the same for up and down parts

```
errorbar(I,V,dV,'x')
ylim([0,6])
```

Yet again we need proper labeling

``` 
title('Dependence of voltage on current')
xlabel('Current (A)')
ylabel('Voltage (V)')
``` 

[errorbar_plot.png]

= Fixing plots font size =

Honestly Matlab is not the best tool to do presentation quality plots, it
requires a quite a bit of messing with their GUI. But at least for a quick font
size fix do the following. This might look like  a magic spell so do read documentation
on **``set``** command.

```
fontSize=24;
set(gca,'FontSize',fontSize );

errorbar(I,V,dV,'x')
title('Dependence of voltage on current')
xlabel('Current (A)')
ylabel('Voltage (V)')
```

[errorbar_with_large_font_plot.png]


= Fitting and data reduction =

Modern experiments generate  a lot of data, while humans  can keep in their
mind only  about 7 things.  So there is no  way to intelligently  work with
large  data set.  We  need  to extract  essential  information  out of  our
data  i.e. do  data  reduction.  Often we  have  some  idea for  analytical
representation/model of  the experimental  data but we  are not  sure about
the  model parameters  values. The  process  of extracting  them is  called
**fitting**.

So let's  look at  above Voltage  vs Current  dependence. According  to the
Ohm's law  such dependence  must be linear  i.e **``V =  R*I``**. Where
**``R``** is constant coefficient called resistance.

So first  we need to  define the model for  Matlab. We **must**  specify an
independent variable name  unless it is called **``x``**, which  is not our
choice  here. Also,  notice  that **``fittype``**  function  wants to  know
independent variable //name// thus **``I``** is surrounded by **``' '``**

```
f=fittype( @(R,I) R*I, 'independent', 'I' )
%______________^ independent variable must be last in the function arguments
%____________^ even though R is supposed to be fixed it is still a variable
%              from fitting point of view
```

Now  our job is simple, we need to specify a guessed value of the R as starting point
for the fit function.

``` param_guessed = [ 4 ]

You might  have a  model with  multiple parameters  than it  is a  bit more
evolved. Please, do read help file about fitting.

Finally, we fit.
```
[fitobject, goodness_of_fit] = fit (I, V, f, 'StartPoint', param_guessed)
```

You will see the following
```
fitobject = 

     General model:
     fitobject(I) = R*I
     Coefficients (with 95% confidence bounds):
       R =       1.291  (0.8872, 1.695)

goodness_of_fit = 

           sse: 2.6340
       rsquare: 0.7050
           dfe: 4
    adjrsquare: 0.7050
          rmse: 0.8115
```

The most interesting part is **``R =       1.291  (0.8872, 1.695)``**, 
where 1.291 is the mean for the extracted parameter ``**R**``
and **``(0.8872, 1.695)``** is the interval within one standard deviation. 
So we get our error bars as well: 
one sigma is equal ``1.695 - 1.291 = 0.4040``

So our conclusion is that **``R = 1.3 +/- 0.4``**.

== Checking the fitting result ==

We should not trust computed fit parameters blindly. For our simple case it
is probably  OK, but  for more evolved  nonlinear fitting  functions result
might  go way  off  if we  specify starting  point  significantly far  from
optimal.

One way  to see check  fit quality  is to look  at Root Mean  Squared error
(**``rmse``**) parameter from above output.  It reflect what is the typical
deviation of the model/fit from experimental points. It should be about the
size  of your  experimental error  bars. If  it is  much smaller  then most
likely you are over fitting (you model has to many free parameters). If it is 
much bigger your fit model is not good enough.


Another way to check consistency of the fit it to look at plotted data and the fit
simultaneously. Let's open yet one more window and plot our data first
```
figure(3)
fontSize=24;
set(gca,'FontSize',fontSize );

errorbar(I,V,dV, 'x')

hold on
```

notice **``hold on``** statement, now next plotting command will draw over presented plot i.e. graphics window //holds on// to already present content.

Now we the plot corresponding to our model/fit. There are multiple way to do it,
the simplest is
```
plot( fitobject )

title('Dependence of voltage on current')
xlabel('Current (A)')
ylabel('Voltage (V)')
```

[fitted_data.png]

This is not a very good fit since the fit line most of the times more than a typical experimental error away.

== More elaborated fit ==

It looks like our  voltage has some offset. We need to  update our model to
include **``Vb``** the biasing voltage term: **``V = R*I +Vb``**. So in our
fit/model should include two unknown parameters **``R``** and **``+Vb``**

```
f=fittype( @(R, Vb, I) R*I+Vb, 'independent', 'I' )
```

We need provide initial
guess for two parameters

```
param_guessed = [ 4, .5 ]
%_________________^           R value
%____________________^^       Vb value
```

The rest of the procedure is the same as above

```
[fitobject, goodness_of_fit] = fit (I, V, f, 'StartPoint', param_guessed)
```

and we see fit parameters **``R``** and **``+Vb``** below
```

fitobject = 

     General model:
     fitobject(I) = R*I+Vb
     Coefficients (with 95% confidence bounds):
       R =      0.9094  (0.5556, 1.263)
       Vb =       1.165  (0.2817, 2.048)

goodness_of_fit = 

           sse: 0.3832
       rsquare: 0.9571
           dfe: 3
    adjrsquare: 0.9428
          rmse: 0.3574
```

Let's see how our new fit looks like

```
figure(4)
fontSize=24;
set(gca,'FontSize',fontSize );

errorbar(I,V,dV, 'x')
hold on

plot( fitobject )

title('Dependence of voltage on current')
xlabel('Current (A)')
ylabel('Voltage (V)')
```

[fitted_data_improved.png]

This seems to be a much better fit.


= Saving you Matlab plots =

Well it nice to have Matlab  plots.  But as soon as you close Matlab they will go away.
How to save them for a future use? Actually, quite easy with **``print``** command. Unlike the name suggest it does not make a carbon copy but actually an electronic one.

For figure saved in png format do
``` print('-dpng', 'V_vs_I.png')

If you find that you figures are to large specify resolution in dots per inch 
with **``-r``** option. For example
``` print('-dpng', '-r50',  'V_vs_I.png')

or if you need a pdf (Matlab usually make very bad pdfs without a lot of tweaking)
``` print('-dpdf', 'V_vs_I.pdf')

The **``-d``** stands for output driver option, i.e png or pdf in above examples.

Make sure that you read documentation for **``print``** command there are some
useful parameters.



= Where to learn more =

At W&M Physics 256 class teaches how scientists use computers. Matlab was a language of
choice for this class. Feel free to read and use
[my Physics 256 materials http://physics.wm.edu/~evmik/classes/2012_fall_practical_computing_for_scientists/]
make sure that you read at least Lecture 02 slides. Also lecture 15 slides discuss
data reduction and fitting.