#                        Celle Ligure, 21/05/2011
#                            Sala consiliare
#                        Matematica in cucina
#   Symposium a cura dell'Associazione culturale-scientifica Cellelab
#
#                   [Ctrl + "+"  per ingrandire]
#
# Premessa  (una barzelletta raccontata dai "fisici")
# Un matematico e un fisico in cucina. Il primo giorno viene loro
# chiesto di cucinare un uovo in camicia. Dispongono di una padella
# attaccata con un filo al soffitto. Il fisico con un coltello stacca
# la padella dal filo, cucina l'uovo, quindi lava la padella e la ripone
# nel mobile accanto. Il matematico con un coltello stacca la padella
# dal filo, ... procendendo allo stesso modo.
# Il secondo giorno viene chiesto loro di cucinare un altro uovo.
# Il fisico prende la padella dal mobile, cucina l'uovo, lava la padella
# e la ripone nel mobile.  Il matematico prende un filo, lega la padella,
# l'appende al soffitto e si rifà al caso precedente.
#
# Spero che, alla fine di questo intervento, avrete idee un po' diverse
# sui matematici ...
#
#                          Carlo Dapueto
#                       In cucina con R
#       http://macosa@dima.unige.it/ingiro/celle/rcucina.htm
#
# R è un programma di statistica, il più usato al mondo (poco in Italia),
# gratuito, che gira su tutte le piattaforme, usabile professionalmente e
# nei vari tipi di scuola, con diversi livelli di impiego, e non solo
# di tipo statistico.
# Vedremo come usarlo per affrontare alcuni problemi da "cucina", ma, nello
# stesso tempo, vorremmo dare un'idea di come lo si può usare in generale.
# Vedremo alla fine come scaricarlo, facilmente, sul proprio computer.
#
# Caratteristica essenziale di R è che tutti i comandi possono essere
# registrati in formato testo, che possono poi essere copiati e riusati
# facilmente, con enormi risparmi di memoria e notevoli facilità d'uso.
# Gli esempi che vedremo sono tutti caricati completamente in questo
# documento (che se vorrete potrete recuperare all'indirizzo sopra
# indicato). Le uscite no!  Se volete aprite il documento ed eseguite
# con dei "copia (dal documento) e incolla (in R)", come farò io,
# i comandi.
# Tutte le righe che incominciano con "#" sono intese come commenti, per
# cui non danno luogo ad alcuna esecuzione da parte del programma.
#
# Partiamo da un cocktail:
## Se mescolo 0.40 litri di succo di arancia e 0.10 litri di liquore a
## 40 gradi, qual è la gradazione alcolica del cocktail ottenuto?
# Il calcolo è semplicissimo (ma se i numeri fossero diversi ...):
(40/100*0.10)/(0.40+0.10)*100
# Ho aperto R, ho copiato le righe precedenti e le ho incollate.
# Ho ottenuto 8:  8% è la gradazione alcolica del cocktail.
# Vediamo un esempio più complesso, rimanendo alle bevande:
## Quanto alcol è contenuto in 200 ml (un bicchiere) di vino dalla
## gradazione del 12% (sappiamo che l'alcol ha densità di 0.79 g/cm³)
12/100*200; 12/100*200*0.79
# Dunque un bicchiere contiene 24 ml di alcol, che corrispondono a 19 g
# Abbiamo visto che posso impostare più calcoli sulla stessa riga;
# basta che siano separati da un punto e virgola
#
# Rimaniamo sui liquidi.  So che 1 litro d'olio pesa 930±10 g.
# Qual è il volume di 1 kg di olio?
m <- (930-10)/1000; M <- (930+10)/1000; v <- c(1/m, 1/M); min(v); max(v)
# Posso assegnare valori a variabili, mediante "<-"; posso mettere in
# una variabile collezioni di valori; posso con i comandi min e max
# trovare il minimo e il massimo di una collezione di valori.
# (sui ricettari e i libri di cucina, o su articoli di tipo medico ...,
#  si trovano modi buffissimi con cui sono riportati i valori numerici) 
#
# Preparo da mangiare. Uso una bilancia.
## Peso un recipiente ottenendo la misura 30 ± 2 g.
## Poi vi verso della farina e con un'altra bilancia, con maggiore
## portata, lo ripeso ottenendo la misura 290 ± 5 g.
## Qual è il peso della farina?
p1 <- 30; e1 <- 2; p2 <- 290; e2 <- 5 
p1a <- p1-e1; p1b <- p1+e1; p2a <- p2-e2; p2b <- p2+e2  
valori <- c(p2a-p1a, p2a-p1b, p2b-p1a, p2b-p1b)
c( min(valori), max(valori) )
# (modificando la terza riga posso operare anche per addizioni,
#  moltiplicazioni, divisioni)
#
# Facciamo un altro semplice esempio, sempre legato ai preparativi:
## Ho 4 botti in cantina, 2 di barbera, 2 di dolcetto. Voglio del
## dolcetto ma non mi ricordo più in quali botti sia. Allora
## assaggio del vino da ogni botte, fino a che trovo quella giusta.
## Qual è il numero medio di assaggi che dovrò fare?
pr1 <- 2/4; pr2 <- 2/4*(2/3) # pr1 prob. trovarlo al 1^ assaggio, pr2 al 2^
1*pr1 + 2*pr2
# Trovo  1.166667 (come è possibile che sia un numero non intero? ...).
# Come posso ottenere il valore in forma frazionaria?
# Si possono caricare librerie con comandi particolari. Una molto
# elementare e molto utile è MASS, che consente di fare, tra
# l'altro, il calcolo frazionario "esatto":
library(MASS)
fractions(pr1 + 2*pr2)
# Troviamo 7/6. L'esempio era elementare, ma si tratta di un comando molto
# utile se si hanno da fare dei calcoli frazionari (come "in cucina").
#
# Un altro esempio, analogo:
## Dal sito web di una rivista di cucina italiana scarico la seguente
## ricetta per le crêpes:    «Per ottenere 25 pezzi occorrono g 250 di
## farina, 4 uova medie, mezzo litro di latte, g 50 di burro e un pizzico
## di sale».    Voglio fare 40 crêpes. Come trasformo la ricetta?
c(250, 4, 1/2, 50) / 25 * 40
# Ottengo  400  6.4  0.8  80  Il problema sono le 6.4 uova, ma ...
# Dal punto di vista del programma quel che è significativo è che
# ho operato sull'intera collezione come se fosse un numero.
#
# Vediamo un altro problema di proporzionalità, ma un po' più complesso:
## «Per preparare della frutta sciroppata ho predisposto 600 g di sciroppo
## (= acqua+zucchero) al 20%. Poi leggo sul ricettario che lo sciroppo deve
## essere al 30%. Quanto zucchero devo aggiungere allo sciroppo che ho già
## preparato?»
#       (è un problema tipico, ma non è un problema "semplice")
# Per risolvere il problema potremmo procedere algebricamente. Vediamo come
# affrontarlo usando il concetto di funzione e tracciando il grafico col
# comando plot:
# ( devo cercare per quale x  (600*20/100 + x)/(600*80/100) = 30/70  dove
#  30/70 è il rapporto tra zucchero ed acqua)
f <- function(x) (600*20/100 + x)/(600*80/100) - 30/70
plot(f,0,200); abline(h=0,v=0, col="blue")
# Dal grafico osservo che la soluzione sta tra 50 e 100 e, con un opportuno
# comando, posso trovarla:
uniroot(f, c(0, 100))$root
#
# Problemi di proporzionalità come questo sono comuni "in cucina"
# (ma sono tutt'altro che facili da affrontare, come risulta da molti
# test di ingresso all'università, anche a Matematica e a Fisica)
# Vediamo una, diversa, questione di proporzionalità:
## Ecco una ricetta per la crema al mascarpone, formalizzata:
#  (per ogni albume metto 2 tuorli, 250 g di mascarpone, …)
Albumi <- function(x) x
Tuorli <- function(x) 2*x
Mascarpone <- function(x) 250*x
Zucchero <- function(x) 80*x
Crema <- function(x) c( Albumi(x),Tuorli(x),Mascarpone(x),Zucchero(x) )
# È una funzione (composta) ad 1 input e 4 output
Crema(1); Crema(5); Crema(7)
#
# Alleggeriamo un po' il discorso con una animazione (senza entrare
# nel merito dei comandi), anche per dare un'idea delle cose che si
# possono fare con R.
# Una "ebollizione":
tic<- function(x) {sec<- proc.time()[3]; while(proc.time()[3]< sec+x) sec<- sec}
view <- c(0,100, 0,70)
col <- function(x) if (runif(1)<1/3) "brown" else { if(runif(1)<0.5) "blue" else "black"}
for (k in 59:0) {tic(0.02); plot(c(view[1],view[2]),c(view[3],view[4]),type="n",xlab="",ylab="")+
abline(h=k+1,col="blue");abline(h=k+1.25,col="blue");abline(h=k+1.5,col="blue")+
if (k>0) for (i in 1:99) for (j in 1:k) if (runif(1)<0.03) points(i,j,col=col(1))}
#     (nota:  ci sono cicli, istruzioni a più righe, ...)
#
# Qui c'è una tabella di composizione alimenti. Vediamo uno dei modi in
# cui possiamo elaborare dati di questo tipo  (consideriamo solo
# proteine, grassi, carboidrati e calcoliamo le calorie corrispondenti):
gr_kcal <- c(4.01,9.03,3.77); etichette <- c("pro","gra","car")
frase <- "Ogni 100 g, g di proteine,grassi,carboidrati e kcal:"
gr<- c(15,4,1); alim<- "Acciughe"
dev.new(width=2.5,height=3); barplot(gr, names.arg=etichette, ylim=c(0,20),main=alim)
paste(alim,frase,sep=". "); gr; round(sum(gr*gr_kcal))
# l'unica cosa che cambio, ogni volta, è  gr<- c(...); alim<- "..."
gr<- c(16,2,0); alim<- "Trippa"
dev.new(width=2.5,height=3); barplot(gr, names.arg=etichette, ylim=c(0,20),main=alim)
paste(alim,frase,sep=". "); gr; round(sum(gr*gr_kcal))
#
gr<- c(8,1,14); alim<- "Piselli"
dev.new(width=2.5,height=3); barplot(gr, names.arg=etichette, ylim=c(0,20),main=alim)
paste(alim,frase,sep=". "); gr; round(sum(gr*gr_kcal))
#
gr<- c(5,5,6); alim<- "Yogurt intero"
dev.new(width=2.5,height=3); barplot(gr, names.arg=etichette, ylim=c(0,20),main=alim)
paste(alim,frase,sep=". "); gr; round(sum(gr*gr_kcal))
#
# Abbiamo visto come dimensionare via via nuove finestre (che poi comunque
# potremmo ridimensionare a mano). Vediamo un diagramma a torta:
dev.new(width=2.5,height=3); pie(gr, labels=etichette, col=c("red","yellow","blue"))
#
# I grafici posso poi (cliccando col pulsante destro sopra essi) salvarli
# come "bitmap", incollarli in Paint, salvarli come GIF, e ...
#
# I dati potremmo averli memorizzati in molti formati e comunque leggerli.
# Ad esempio se li avessi memorizzati nel file "alimenti.txt" collocato
# in qualche cartella, col comando "Change Dir" del menu "File" potrei
# spostarmi su tale cartella e mettere i dati in una variabile:
dati <- read.delim("alimenti.txt",header=FALSE,sep =",")
dati
# Se non ho il file precedente sul computer (se non l'ho salvato o se
# non l'ho nella cartella dei file temporanei) posso mettere l'indirizzo
# completo:
indirizzo <- "http://macosa.dima.unige.it/ingiro/celle/alimenti.txt"
dati <- read.delim(indirizzo,header=FALSE,sep =",")
# Posso poi elaborali in vario modo. Vediamo, solo, come elaborare le proteine:
min(dati$V2); max(dati$V2); sort(dati$V2); mean(dati$V2); median(dati$V2)
stem(dati$V2)
# L'ultima rappresentazione "stem-and-leaf" (diagramma a "gambo e foglia")
# è il più usato strumento di rappresentazione grafica di dati statistici
#
# 0 | 0224     È un istogramma realizzato in formato testo che
# 0 | 588      mi dà anche le informazioni numeriche sui singoli
# 1 |          dati elaborati. In questo caso ho  0,2,2,4 tra 0 e 5 escluso,
# 1 | 5668     5,8,8 tra 5 e 10 escluso,  nessun dato tra 10 e 15 escluso,
# 2 | 00       15,15,16,18 tra 15 e 20 escluso,  20,20 tra 20 e 25 escluso.
#
# Vediamo un'altra situazione, in cui ho una grandezza legata ad un'altra
# mediante una funzione di cui non conosco o non ricordo l'espressione: qui
# (temperatura di ebollizione dell'acqua in funzione della pressione atm.)
# traccio i punti
t <- c(0,  25,  40,50,  60, 65, 70, 75, 80, 85, 90, 95, 100,105,110)
p <- c(4.6,23.8,54,92.5,156,194,238,290,355,434,526,634,760,906,1075)
plot(p,t); abline(h=0,v=0,col="blue")
# li congiungo con dei segmenti
lines(p,t)
# 100 gradi corrispondono a 760 mmHg
abline(h=100,lty=3); abline(v=760,lty=3)
# troviamo graficamente la altezza a cui corrisponde la pressione 720
abline(v=720,lty=3)
# usiamo il comando SPLINE che approssima i punti con una "buona" curva
spline(p,t,xmin=700,xmax=740)
# vediamo che a 720 e' associato 98.49183
# se ho difficoltà di individuazione basta che batta:
spline(p,t,xmin=720,xmax=720)
# facciamo uno zoom e vediamo meglio la cosa
plot(p,t,xlim=c(600,800),ylim=c(90,100)); abline(h=0,v=0,col="blue")
lines(p,t); abline(v=720,lty=3)
abline(h=98.49,lty=3)
# Il segmento, giustamente, sta un po' sotto al punto stimato con la spline
# Abbiamo illustrato con un esempio, ma sono molti gli aspetti (anche
# in ambito "culinario") che possono essere stimati con tecniche di questo
# tipo.
#
# Un altro esempio "in cucina":  Come si fa la pasta al dente a 3000 m?
# Devo tener conto dei cambiamenti di pressione, che, medianemente,
# sono indicati dalla allegata tabella
# So che:  760 mmHg = 1013.25 hPa
700*760/1013.25
# Ottengo 525 hPa; dal grafico precedente ottengo che bolle a 90°:
# devo cuocere la pasta in una (opportuna) pentola a pressione!
#
# Vediamo, rapidamente, due esempi legati all'analisi matematica.
# Anche qui si potrebbero fare vari spot legati a situazioni abbastanza
# standard, come il tempo che impiega un termometro da forno estratto
# da un forno caldo a raffreddarsi (si fanno di leggere i valori, li si
# rappresentano graficamente e si "trova" che l'andamento è esponenziale),
# o, analogamente, il tempo che impiega un alimento messo in un forno
# a raggiungere la temperatura fissata.
# Ma si tratta di ragionamenti che comporterebbero via un po' di tempo.
# Vediamo, dunque, due semplici esempi, sempre, in qualche modo, legati
# alla "cucina":
#
## «Riempo con acqua corrente un recipiente non cilindrico. In particolare
## un contenitore conico di altezza 2 e raggio 1. Qual è la velocità di
## variazione di h rispetto a V?»     So che h = (12*V/π)^(1/3):
D(expression((12*V/pi)^(1/3)),"V")
# Abbiamo fatto la derivata.
#
## «Negli anni intorno al 2000 la nazione con maggiore consumo individuale
## di carne di maiale è stata la Danimarca. Sappiamo che nel 1994 il consumo
## procapite è stato di 66.7 kg con una deviazione standard di 28.1 kg.
## Se avessimo preso, a caso, 25 danesi, qual è la probabilità che la media
## annua del loro consumo di carne di maiale avesse superato i 68 kg?»
# Senza entrare in dettagli, osserviamo che il loro consumo medio (essendo 25
# sufficientemente grande) per il "teorma limite centrale" ha andamento circa
# "normale" (ossia "gaussiano"), di media 66.7 e deviazione standard 28.1/√25
m <- 66.7; s <- 28.1/sqrt(25);  d <- function(x) dnorm(x,mean=m,sd=s)
plot(d,m-3*s,m+3*s); abline(h=0,col="blue"); abline(v=c(m,m-s,m+s),lty=3)
integrate(d, 68, Inf)
# Si è visto come calcolare un integrale definito in un intervallo infinito
#
# Finiamo, sperando di aver "toccato" gli aspetti principali di come la
# "matematica" può entrare "in cucina". Chi è interessato a sperimentare
# l'uso di R può andare qui e vedere come scaricarlo (ed usarlo).
#  (chi fosse interessato ad esercizi che intrecciano matematica e scienze
#   biologiche clicchi qui).
#
# Ciao ...
#
# -------------
# [QUI chi vuole, dopo aver provato (o prima di farlo) può vedere le uscite]