#
#   - - - - - - - -   INCONTRO N. 1 - Primi "comandi"   - - - - - - - -
#
# Vediamo, in dettaglio, la spiegazione di alcuni comandi e alcuni esempi
# di applicazione.
#
# Posso calcolare termini qualunque, rispettando le priorita' tra le
# operazioni. Un ";" consente di separare due comandi sulla stessa riga  
2+4*4/sqrt(3)^2; (2+4)*4/sqrt(3)^2; 2+4*(4/sqrt(3))^2; 2+(4*4/sqrt(3))^2
#
# Posso mettere dei valori in variabili, ossia in espressioni inizianti
# con una lettera; per assegnare un valore ad una variabile occorre usare
# una freccia, ossia <- o ->. 
# Posso assegnare lo stesso valore contemporaneamente a più variabili,
# concatenando più frecce.
# A differenza di altri programmi, per R nomi che differiscono per le
# dimensioni dei caratteri sono interpretati come differenti.
a <- 7; 5 -> k; a+k; a*k; c <- pippo <- a*k+a*k+1/2; pippo; c+1; c+K
#
# Vengono visualizzate 7 cifre ma i risultati dei calcoli sono memorizzati
# con (circa) 16 cifre buone; vediamo cosa accade a 3 numeri, collezionati
# in un'unica variabile:
x <- c(10, 45, 5/10^5); x/3; print(x/3, digits = 16)
#
# Occorre tener presente che nel caso del risultato di un calcolo possono
# intervenire errori di arrotondamento ed aversi meno cifre "buone":
x <- 5555.1251; y <- 5555; x-y; print(x-y, digits = 16)
#
# Schiacciando i tasti freccia [^] e [v] si possono rivedere i precedenti
# comandi: cosa COMODA per riusare, correggere o modificare comandi già
# usati. Prova a farlo.
#
# Se sei in rete (o se ne hai fatto un uso recente da rete) puoi guardare
# uno dei vari help; ad es. selezionare "R functions (text)", mettere "sin"
# come parola da cercare; si apre "Trig {base}" (da qui, volendo, seleziona
# ulteriormente "Math").  Puoi usare anche "Search Engine & Keywords" attiva-
# bile da "Html Help" e, soprattutto, "Base" e "Graphics", in "Packages",
# sempre attivabile da "Html Help"; si tratta, comunque, di help non facili
# per il principiante. Battendo il comando  help(sin)  si entra direttamente
# nell'help di sin.
# Non da rete, puoi azionare "Search Help"; se poi batti "sin" hai l'elenco
# delle voci dell'help in cui si parla di sin (tra cui Trig).
# L'help online che stiamo realizzando vuole mettere a disposizione un "aiuto"
# semplice, usabile a scuola.
#
# Per la radice cubica, quinta,..., n-ma per n dispari si può usare ^(1/n)
# che però è definita solo per input positivi; altrimenti si usano le
# funzioni sign (segno) e abs (valore assoluto); es.:
x <- -8; x^(1/3); sign(x)*abs(x)^(1/3)
# NaN sta per "non è un numero"
#
3+sqrt(0.25); 3+sqrt(-0.25)
# Il programma non calcola la radice quadrata di numeri negativi (e non
# fornisce quindi gli strani risultati che danno Derive e Maple, per i quali
# x*√x/√x equivale ad x anche per x negativo:
# vedi   http://macosa.dima.unige.it/om/ind/comples16.htm
# Opera anche sui numeri complessi solo se compare una "i":
3+sqrt(-0.25+0i)
#
# Ovviamente, non tutto puo' essere calcolato:
factorial(14); factorial(170); print(factorial(170),16); factorial(171)
# Come si è visto, delle uscite intere sono stampate anche più di 7 cifre
#
# Altro esempio di calcolo di funzione ad output interi: il coeff. binomiale
choose(42,21); choose(54,27); print(choose(54,27),16)
# Ma per coeff. binomiale con più di 16 cifre si possono ottenere valori non
# attendibili. Ad es. il coeff. binom. (76,38) è 6892620648693261354600 ma
# c.bin.(76,38) fornisce  6.892621e+21, 6892620648693271363644
choose(76,38); print(choose(76,38),22)
#
# Vediamo come trovare i divisori di un numero naturale; e vediamo come usare
# un ciclo for e un'istruzione if (nella condizione occorre usare == invece
# di =); ricordiamo che floor è la parte intera.
n <- 68; i<-0; for(k in 1:n) if (n/k==floor(n/k)){i <- i+1;print(c(i,"° divisore",k))}
#
# Miglioriamo l'aspetto estetico delle uscite: con l'opzione "quote=FALSE"
# si evita la stampa delle "virgolette"; con "paste" si incollano le uscite
# separandole con il separatore indicato in "sep" (qui con "", ossia niente).
# Scriviamo poi l'istruzione su due righe: il programma, automaticamente,
# interpreta come un'unica istruzione quella spezzata su più righe.
n <- 68; i<-0; for(k in 1:n) if (n/k==floor(n/k)) {i <- i+1; 
        print(c(paste(i,"° divisore",sep=""),k),quote=FALSE)}
#
# Vediamo un altro esempio di funzione:
# Un piccolo gelataio produce e vende un solo tipo di cono, con i seguenti
# costi: 0.50 euro di spese incoporate (per latte, cacao, cialde, ...) e
# 2000 euro di spese fisse mensili (per locali, energia elettrica, ...)
SpeseFisse <- 2000;  SpeseIncorporate <- 0.50
# Quanto gli costa un gelato se al mese ne vende N?
CostoUnitario <- function(N) SpeseIncorporate + SpeseFisse/N
# Se ne traccio il grafico per N tra 0 e 5000 ottengo:
plot(CostoUnitario,0,5000)
#
# Il programma sceglie da solo la scala verticale. Usando ylim=c(.,.)
# posso scegliere dove far variare le coordinate verticali; suppongo che
# mi interessi vedere quando il costo unitario è inferiore a 10 euro
plot(CostoUnitario, 0,5000, ylim = c(0,10)); abline(h=0, v=0)
#
# Nel grafico precedente c'è anche, tratteggiata, la retta orizzontale
# passante per 2. Se il prezzo di un gelato è 2 euro, essa
# rappresenta tale prezzo. È stata tracciata con:
abline(h=2,lty=3)
# È stata tracciata anche la retta verticale che passa per il punto del
# grafico con ordinata 2. Ecco come è stata trovata la ascissa:
# da 2 = SpeseIncorporate+SpeseFisse/n    ricavo
n <- SpeseFisse/(2-SpeseIncorporate); n
# e come è stata tracciata la retta verticale di tale ascissa:
abline(v=1333.333,lty=3)
#
# In questo caso era facile risolvere l'equazione. In generale
# se l'equazione corrisponde al punto in cui il grafico di una
# funzione continua attraversa una retta si può trovare la
# soluzione con un metodo numerico (l'idea è il semplice metodo
# di bisezione, che poi vedremo come implementare direttamente).
dove <- function(n) CostoUnitario(n)-2
plot(dove,0,5000,add=TRUE,col="red")
uniroot(dove, c(1000, 2000), tol = 1e-10)
# Vengono stampati tutti i dettagli del calcolo. Se ci interessa solo
# la soluzione, specifichiamolo con l'aggiunta di $root:
uniroot(dove, c(1000, 2000), tol = 1e-10)$root
# Sopra abbiamo anche "aggiunto" il grafico della funzione di cui
# abbiamo trovato lo "zero".
#
# Vediamo altre cose utili.
# Ecco come conviene salvare un'immagine. Si clicca su di essa e, dal menu
# che si apre cliccando, la si salva come BitMap e la si incolla su Paint
# (o un altro programma di grafica), si aggiunge quello che si vuole e la
# si ridimensiona. A questo punto o la si copia e salva nel documento in
# cui si sta operando o, da Paint, la si salva in formato GIF, che occupa
# poco spazio e consente di operare successive modifiche, anche dettagliate,
# con Paint (nelle ultime versioni di Windows, conviene, prima di salvarlo,
# farsi una copia del file, quindi memorizzarlo come GIF, attendere un
# messaggio che dice che l'immagine non è salvata perfettamente, quindi
# incollare l'immagine memorizzata e fare Save: la immagine a questo punto
# è perfetta). Provate a salvare l'ultimo grafico.
#
# Per visualizzare le coordinate di punti su un grafico dò il comando
# locator(N) con N = numero dei punti e faccio dei clic su di essi, e
# ottengo approssimazioni delle loro ascisse e delle loro ordinate. Con
# locator(1)$x ho la visulizzazione solo del valore della ascissa, con
# locator(1)$y solo quello dell'ordinata). Proviamo ad usare locator.
#
# Se voglio (ad es. per conservare la vista di grafici precedenti), posso
# aprire una o più nuove finestre grafiche usando il comando dev.new()
# (nuovo dispositivo)
# Un esempio, al di là della sua significatività:
plot(CostoUnitario,0,5000); abline(v=0,h=0,col="red")
dev.new()
plot(CostoUnitario,-1000,1000); abline(v=0,h=0,col="blue")
#
# Vedremo le prossime volte come memorizzare/caricare file di comandi,
# file di dati nei più svariati formati, ..., come realizzare
# grafici tridimensionali, ...
#
# Vediamo solo uno dei vari modi in cui i dati possono essere recuperati dal
# computer o, in rete, senza dei copia/incolla; ad esempio col comando:
readLines("t_arrivi.txt", n=5)
# In prima battuta dovrebbe comparirci un messaggio in cui siamo avvisati
# che il file non viene trovato. Il file è nella stessa cartella in
# cui è il file che stiamo leggendo. Basta che utilizziamo il comando
# Change Dir dal menu Edit per spostare la cartella di riferimento (posso
# anche copiare il percorso dalla finestra del browser e incollarlo). Se
# poi ribattiamo il comando, esso viene caricato:
readLines("t_arrivi.txt", n=5)
# Se siamo in rete possiamo mettere direttamente il percorso per arrivare al
# file:
readLines("http://macosa.dima.unige.it/om/prg/R/agg/t_arrivi.txt", n=5)
# Leggiamo le prime 5 righe del file. Vediamole tutte (in questo caso
# sono poche):
readLines("t_arrivi.txt")
# Col seguente comando metto in "DATI" il file saltando con skip 1 riga
DATI <- scan("t_arrivi.txt", skip=1)
DATI
# Sono gli stessi dati caricati all'inizio dell'incontro:
summary(DATI)
#
# Facciamo due ultimi esempi:
#
# Calcoliamo i valori per le durate di 5, 10, 15, …, 50 anni
# della crescita di 100 € secondo i tassi annuali di interesse composto,
# capitalizzato mensilmente, pari al 6%, all'8%, al 10% e al 12%, e
# rappresentiamoli graficamente.
i <- c(6,8,10,12)/12/100; p <- function(x) round(x,2)
for (a in (0:10)*5) { m <- a*12; print( p (c(a,m,100*(1+i)^m))) }
# i grafici:
plot(c(0,50),c(0,10e3),type="n",xlab="", ylab="")
abline(h=0,col="brown"); abline(v=0,col="brown")
abline(h=seq(0,1e4,1e3),v=seq(0,50,10),lty=3,col="grey50")
axis(1,pos=0,label=FALSE,col="brown"); axis(2,pos=0,label=FALSE,col="brown")
for (j in i) for (a in (0:10)*5) {m <- a*12; points(a,100*(1+j)^m)}
f <- function(x) 100*(1+j)^(x*12)
for (j in i) curve(f,add=TRUE)
#
# Ultimo esempio. Una animazione (usa una finestra non troppo piccola):
yy <- c(-1,-0.75,-0.5,-0.25,0,0.25,0.5,0.75,1); xx <- c(0,0.25,0.5,0.75,1)
ce <- function(x) sqrt(1-(x^2))
n <- 1000; for (i in 0:n)                                             #inizio
{{plot(c(-1,1),c(0,1),type="n",xlab="",ylab="",asp=1)}
{abline(h=xx,col="grey");abline(v=yy,col="grey");abline(h=0,col="blue")}
{curve(ce, add=TRUE,col="blue");x <- cos(pi/n*i);y <- sin(pi/n*i)}
{x1 <- c(-1,x,1); y1 <- c(0,y,0);lines(x1,y1);points(x*0.9,y*0.9)}}   #fine
# Non entriamo nei dettagli dei comandi. Questo esempio serve solo a dare una
# idea delle cose che si potrebbero approfondire.
#
# Torniamo indietro (<- clicca).