torstai 31. maaliskuuta 2016

Machine learning: Bagging ja puumallit

Esittelen koneoppimiseen kuuluvan bagging-menetelmän, jolla puumallien ennustusten varianssia pystytään vähentämään. Puumallien yleinen ongelma on ylisovittaminen (Overfitting), jonka takia mallit ennustavat hyvin harjoitusdataan, mutta huonosti uuteen dataan. Ne siis löytävät harjoitusdatasta rakenteita, jotka eivät kuitenkaan ole yleistettävissä muuhun dataan. Alla kuvassa malli on ylisovittanut ja hakenut punaisia pisteitä sinisten joukosta. Musta viiva antaa kuitenkin parhaan ennustustuloksen, kun mallia kokeillaan uuteen dataan. Punaiset pisteet sinisten joukossa eivät siis tässä muodossa toistu usein, kun otetaan uusia otoksia perusjoukosta.

Lähde: https://kaggle2.blob.core.windows.net/competitions/kaggle/2489/logos/front_page.png
Puumalli on mahdollista sovittaa täydellisesti harjoitusdataan, jolloin malli ennustaisi virheettömästi harjoitusdatan. Tämä johtuu siitä, että puumalli voi aina tehdä uuden oksan, kunnes kussakin oksan lehdessä on enää 1 havainto jäljellä. Esim. luottoriskiluokituksiin liittyen tämä tarkoittaisi sitä, että jokainen yksilö luokiteltaisiin loppuun asti, eli kukin yksilö kulkisi puussa omaan lehteensä, jossa ei olisi ketään muuta. Esimerkiksi kaksi lainanhakijaa, joista toisella olisi todettu maksukyvyttömyys, voitaisiin puumallilla jakaa eri lehtiin, vaikka heidän ainoa eronsa olisi yksi euro säästötilillä. Tällöin puumalli veisi hakijat samaa reittiä viimeiseen oksaan, jossa se tekisi jaon säästötilin rahasumman mukaan. Jos tällaista puuta käyttäisi luotonantoon, niin näistä kahdesta käytännössä identtisestä hakijasta vain toinen saisi lainan. Ylisovittamisen estämiseksi puumalleihin voi asettaa erilaisia ehtoja, kuten havaintojen minimimäärät lehdissä ja oksien maksimimäärät mallissa. Lisäksi ns. boostrap aggregating (bagging) on tehokas keino välttää ylisovittaminen.

Bagging

Sen sijaan, että käyttäisimme yhden puumallin harjoitusdatasta muodostamaa ennustetta, voimme ottaa harjoitusdatasta monta otosta ja muodostaa niiden pohjalta useita puumalleja. Kun näillä puumalleilla sitten ennustetaan testidataa, niin käytämmekin yksittäisten mallien sijaan niiden ennusteiden keskiarvoa. Tällä tavoin saamme luotua vakaamman ennusteen, jolloin yksittäisten mallien ylisovittuminen ei ole enää yhtä suuri ongelma (Baggingin avulla pienet muutokset datassa eivät aiheuta suuria muutoksia ennusteissa). Havainnollistan baggingin vaikutusta simuloimalla datan, jossa on kaksi muuttujaa y ja x. Lisään y:n vähän "kohinaa", jotta näemme, kuinka lähelle mallin ennusteet menevät, kun kohina häiritsee puumalleja. Muodostan sitten datasta harjoitusdatan (60% havainnoista) sekä testidatan (40% havainnoista). Alla on kuva testsetistä, johon ennustukset tehdään. Kuvassa oleva sininen viiva kuvaa parasta mahdollista ennustetta, koska se on simuloidun datan funktio ilman hajontaa aiheuttavaa satunnaista kohinaa.
R-koodi:

x <- sort(runif(300, 0, 10))
y <- sin(2.5*x)/x + rnorm(300,0,sd=0.10) 
data <- data.frame(y,x)
t <- sort(sample(nrow(data),0.6*nrow(data)))
trainingset <- dataa[t,]
testset <- dataa[-t,]
plot(testset$x, testset$y)
realvalue <- sin(2.5*x)/x
lines(x, realvalue,col="blue",lwd=2)
legend("topright",c("sin(2.5*x)/x",
"sin(2.5*x)/x + rnorm(300,0,sd=0.10)"),col=c("blue","black"),lty=c(1,NA),pch=c(NA,1),lwd=c(2,1))

Seuraavaksi muodostetaan puumalli käyttäen koko trainingset-dataa ja lasketaan sen antamien ennusteiden perusteella RMSE (Root Mean Squared Error, lisää infoa: RMSE). RMSE:lla verrataan ennusteiden osuvuutta oikeisiin arvoihin. Mitä pienempi arvo, sitä paremmin malli on ennustuksissaan onnistunut.
R-koodi:

library(rpart)
m0 <- rpart(y~.,data=trainingset,control=rpart.control(minsplit=5))
testset$pre0 <- predict(m0,testset) 
RMSE <- sqrt(mean((testset$y - testset$pre0)^2))
print(RMSE)
plot(testset$x, testset$y) # Uusi kuva, jossa m0:n ennusteet punaisella
realvalue <- sin(2.5*x)/x
lines(x, realvalue,col="blue",lwd=2)
lines(testset$x,testset$m1,col="red")
legend("topright",c("sin(2.5*x)/x",
"sin(2.5*x)/x + rnorm(300,0,sd=0.10)","m0 (RMSE=0.151)"),col=c("blue","black","red"),lty=c(1,NA,1),pch=c(NA,1,NA),lwd=c(2,1,1))




Ensimmäisen puumallin RMSE arvo on 0.151. Seuraavaksi otan traininigset-datasta viisi 20%:n satunnaisotosta (palautuksella), joista jokaisella muodostan oman puumallin. Ennustan näillä malleilla sitten testsettiin ja lasken RMSE:t. Sitten lasken mallien ennusteista bagging-menetelmän mukaisesti keskiarvot ja lasken uuden RMSE:n.
R-koodi:

a <- sample(nrow(trainingset),0.2*nrow(trainingset))
b <- sample(nrow(trainingset),0.2*nrow(trainingset))
c <- sample(nrow(trainingset),0.2*nrow(trainingset))
d <- sample(nrow(trainingset),0.2*nrow(trainingset))
e <- sample(nrow(trainingset),0.2*nrow(trainingset))
data1 <- trainingset[a,]
data2 <- trainingset[b,]
data3 <- trainingset[c,]
data4 <- trainingset[d,]
data5 <- trainingset[e,]
m1 <- rpart(y~., data=data1, control=rpart.control(minsplit=5))
m2 <- rpart(y~., data=data2, control=rpart.control(minsplit=5))
m3 <- rpart(y~., data=data3, control=rpart.control(minsplit=5))
m4 <- rpart(y~., data=data4, control=rpart.control(minsplit=5))
m5 <- rpart(y~., data=data5, control=rpart.control(minsplit=5))
pre1 <- predict(m1,testset)
pre2 <- predict(m2,testset)
pre3 <- predict(m3,testset)
pre4 <- predict(m4,testset)
pre5 <- predict(m5,testset)
RMSE1 <- sqrt(mean((testset$y - pre1)^2)) #lasketaan kullekin mallille RMSE-arvot
print(RMSE1)

Yllä olevassa kuvassa näkyvät m1-m5:n ennusteet sekä niiden yhteiset mAvg:n ennusteet. Pienimmän RMSE:n antaa malli mAvg, eli bagging on parantanut ennusteita verrattuna yksittäisiin malleihin m1-m5. Alla on vielä yllä oleva kuva ilman sinistä viivaa, jotta bagging:n merkitys näkyy selkeämmin. Kuvasta näkee, kuinka vaihtelevia yksittäisten puumallien ennustukset ovat, vaikka ne ovat kaikki peräisin saman trainingsetin arvoista. Baggingin avulla näistäkin ailahtelevista puumalleista saadaan luotua vakaa ennustus (ennusteet ovat vakaita, koska baggingin takia tapahtuu konvergoitumista, eli lisäämällä otosten määrää trainingsetistä mAvg:n tulokset eivät radikaalisti poikkea nykyisestä), koska keskiarvotus jättää yksittäisten puumallien rajut ylisovittamiset vähemmälle huomiolle. Lisäksi, jos enemmistö puumalleista on ennusteessaan melko oikeassa, niin nämä ennusteet saavat enemmän painoarvoa (esim. katso kuvaa, kun x=3). Tietysti, jos enemmistö puumalleista olisi väärässä, niin myös huonot ennusteet painottuisivat. Kuitenkin Ensemble learning -menetelmät perustuvat ajatukseen, että keskinkertaisista malleista saadaan luotua hyviä yhdistämällä, ja pahimmassakin tapauksessa yhdistetty malli on yksittäisen mallin veroinen, eli ei muutu sitä huonommaksi. Kahden mallin yhdistäminen voi olla hyvä ajatus erityisesti, jos kummankin mallin ennusteet poikkeavat toisistaan ja mallit ovat melko varmoja oikeissa ennusteissaan, kun taas epävarmoja ollessaan väärässä. Tällaiset mallit yhdistämällä niiden oikeaksi menneet ennusteet pysyvät sellaisina myös keskiarvotusten jälkeen ja lisäksi osa toisen mallin epäonnistumisista korjaantuu, koska toinen malli ennusti ne oikein. Näin väärien ennustusten määrä onkin yhdistetyssä mallissa pienempi kuin yksittäisissä. Lisää tieto löytyy mm. täältä

"A necessary and sufficient condition for an ensemble of classifiers to be more accurate than any of its individual members is if the classifiers are accurate and diverse. An accurate classifier is one that has an error rate of better than random guessing on new x values Two classifiers are diverse if they make different errors on new data points"

Lähde: http://web.engr.oregonstate.edu/~tgd/publications/mcs-ensembles.pdf

Yhdistetty malli mAvg on tällä hetkellä m0:n veroinen ennustuksissaan, joten baggingista ei ole ollut käytännössä vielä huomattavaa hyötyä. Seuraavaksi otan viisi uutta satunnaisotosta trainingset-datasta ja muodostan uudet puumallit, joiden ennusteet yhdistän m1-m5:n kanssa. R-koodi on vastaava kuin aiempi. Nyt RMSE on 0.141. Alla kuva, jossa vanha ja uusi malli.
Alla vielä kuva koneoppimiseen liittyvän Random forest -menetelmän ennusteista. RMSE on tälle mallille 0.119. Random forest -menetelmässä algoritmi muodostaa useita puumalleja (tässä 5000 kpl) ja keskiarvottaa niiden ennusteet (luokittelevassa ennustamisessa se käyttää ennusteiden moodia). Puumallit tehdään bagging-menetelmällä, mutta myös muuttujista käytetään vain osaa kussakin puumallissa (muuttujat kuhunkin jakoon (split) valitaan sattumanvaraisesti). Oletuksena R:ssa muuttujia "n" käytetään kussakin jaossa neliöjuuri(n) (luokittelu) tai n/3 (regressio) määrä. Random forest on kätevä menetelmä, jolla saa parempia ennusteita kuin yksittäisillä puumalleilla. Random forest ei kuitenkaan juuri paranna ennusteita verrattuna logistiseen regressioon, kun ennustetaan binääristä selitettävää muuttujaa. Random forest löytää kuitenkin muuttujien välisiä yhteisvaikutuksia (Interactions) ja on siten avuksi logistisen regression rinnalla. Näiden mallien yhdistämisestä kirjoittelen jatkossa.




maanantai 28. maaliskuuta 2016

Lineaarinen regressio ja Interaction-termit

 Lineaarinen regressio sopii tilanteisiin, joissa haluamme tietää, kuinka eri muuttujat vaikuttavat selitettävään muuttujaan määrällisesti. Lineaarisella regressiolla voi esimerkiksi tutkia, miten ihmisten tulotaso vaikuttaa heidän ostamansa auton hintaan. Lineaarisella regressiolla voi laskea myös todennäköisyyksiä jollekin binääriselle tapahtumalle, mutta logistinen regressio sopii tällaiseen huomattavasti paremmin mm. koska logistinen regressio pakottaa todennäköisyyden 0:n ja 1:n välille, kun taas lineaarinen regressio ennustaa myös ulkopuolen (< 0 & >1) arvoja. Mikäli ratkaisee tämän ongelman asettamalla alle nollan arvot nollaksi ja yli yhden arvot yhdeksi, niin syntyy äärimmäisiä 0:n ja 1:n todennäköisyyksiä, eli malli ennustaisi tapahtumia täysin varmoiksi. Lineaarinen regressiomalli rikkoo myös esim. virhetermien normaalisuusehdon, joten logistinen malli on parempi vaihtoehto. Logistisen regressiomallin ennusteet  käyttäytyvät pehmeämmin kuin lineaarisen regressiomallin ("ei lineaarisesti"), kun lähestytään todennäköisyyksien ääripäitä. Alla olevasta kuvasta näkyy, kuinka logistisen mallin todennäköisyydet pehmenevät ääripäissä.

Lähde: http://www.appstate.edu/~whiteheadjc/service/logit/intro.htm

Vaikka lineaarinen regressio ei toimi logistisen mallin veroisesti todennäköisyyksien mallintamiseen, on se hyödyllinen muihin tarkoituksiin. Sen avulla voi saada intuitiivista tietoa nopeasti eri muuttujien välisistä mahdollisista yhteyksistä. Lineaarista mallia käytettäessä on tärkeää, että selittävän ja selitettävän muuttujan välillä vallitsee lineaarinen yhteys. Mikäli yhteys ei ole lineaarinen, on mahdollista käyttää esim. neliöjuuri-muunnettua muuttujaa normaalin muuttujan kanssa, jolloin epälineaarisia yhteyksiä voidaan mallintaa. Tällaisessa tapauksessa täytyy kuitenkin muistaa tarkastella, onko tämä teoreettisesti oikein. Esim. tulotason kasvu voi vaikuttaa ostetun auton hintaan nostavasti tiettyyn pisteeseen "x" saakka, kunnes kasvu hiipuu tai jopa laskee. Tällöin suoran viivan vetäminen ei kuvaa yhteyttä oikein, koska mallin mukaan auton hinta kasvaa vielä, kun menemme kriittisen tulotason "x" yli. Tällaisessa tapauksessa voi harkita selittävän muuttujan luokittelua, jos haluaa pitää mallin helpommin tulkittavana, eli ei lisäile epälineaarisia yhteyksiä löytäviä muuttuja-muunnoksia. Yksi regressiomallien hienouksia on, että voimme tutkia, kuinka muuttuja "x vaikuttaa selittettävään muuttujaan "y" kun kontroillaan muuttujaa "z". Eli voimme esim. tutkia, kuinka tulotaso vaikuttaa auton hintaan, kun pidämme sukupuolen vakiona. Kontrollointi tapahtuu lineaarisessa mallissa seuraavasti:

Haluamme tietää mallin (y=x+z) kertoimen x:lle. Ensin ratkaistaan x:n vaikutus y:hyn regressoimalla x z:lla ja tallentamalla residuaalit. Tämän jälkeen y regressoidaan näillä residuaaleilla. Poistamme siis ensin x:sta z:n vaikutuksen ja sitten vasta regressoimme y:n tämän x:n, mutta z:sta riippumattoman osan kanssa. Vastaavalla tavalla saamme myös z:n kertoimen. Yksi regressiomallin yleisistä ongelmista (Omitted variable bias) aiheutuu juurikin tämän kontrolloinnin puutteesta, eli jos tärkeä muuttuja jätetään mallista pois, niin muuttujien regressiokertoimet voivat olla vääristyneitä. Esim. x:n kerroin voi olla suurempi kuin mitä se olisi jos olisimme ottaneet z:n vaikutuksen huomioon.

Alla vielä vastaava esimerkki R-koodilla:

sukupuoli <- (rbinom(100,1,0.5)) #luodaan sukupuoli-muuttuja Bernoullin jakaumasta
tulot <- rnorm(100,15000,100)+150*sukupuoli # luodaan toinen muuttuja, joka korreloi sukupuolen kanssa
y=10000+(0.1*tulot)+(500*sukupuoli) # y kuvaa auton hintaa. Tästä näemme oikeat regressiokertoimet.

#Ratkaistaan ensin tulojen regressiokerroin

regTs <- lm(tulot~sukupuoli)
residuaalitT <- residuals(regTs)
m1 <- lm(y~residuaalitT)
round(coef(m1),1) # Saadaan kerroin 0.1. Round-komento pyöristää tässä luvut yhteen desimaaliin.
regSt <- lm(sukupuoli~tulot) # Ratkaistaan vastaavasti sukupuolen regressiokerroin 
residuaalitS <- residuals(regSt)
m2 <- lm(y~residuaalitS)
round(coef(m2),2)# Saadaan 500

sum(residuaalitT*y)/(sum(residuaalitT^2)) #Regressiokertoimet saadaan laskettua myös näin
sum(residuaalitS*y)/(sum(residuaalitS^2)) 

sum(y-0.1*tulot-500*sukupuoli)/100 # Intercept saadaan lopuksi tästä kaavasta
#Lopputulos on siis y=10000+(0.1*tulot)+(500*sukupuoli)

Omitted variable bias tulee ilmi, kun jätämme sukupuolen pois tarkastelusta:

m3 <- lm(y~tulot)
round(coef(m3),1)
#saadaan tulojen kertoimeksi 0.3

Koska tulot ja sukupuoli korreloivat positiivisesti (cor=0.63) ja sukupuoli vaikuttaa auton hintaan positiivisesti, niin tulojen kerroin on liian suuri, jos jätämme sukupuolen pois mallista. Vastaavasti, jos sukupuoli vaikuttaisi auton hintaan negatiivisesti (y=10000+(0.1*tulot)-(500*sukupuoli)), niin tulojen kerroin olisi alaspäin vääristynyt -0.1.

Interaction-termit

Yllä olevassa esimerkissä tulojen kasvu vaikuttaa auton hintaan kullakin sukupuolella samoin . On kuitenkin mahdollista, että eri sukupuolten väliset preferenssit eroavat ja toinen sukupuoli käyttää tulojen noususta suhteellisesti enemmän  rahaa autoonsa. Tällainen yhteys voidaan mallintaa lisäämällä malliin ns. Interaction-muuttuja, jolloin tulojen muutos voi vaikuttaa eri tavoin auton hintaan eri sukupuolilla. Interaction-muuttuja lisätään R:lla seuraavasti: 

m4 <- lm(y~tulot+sukupuoli+tulot:sukupuoli)

Nyt tulojen kerroin riippuu sukupuolesta. Esim. sukupuolella "1" tulojen kerroin voi olla 0.5 ja sukupuolella "0" 0.7. 

Alla on kuva datasta, jossa tulot vaikuttavat ostetun auton hintaan (y) eri tavoin. Selvästi sukupuoli "0" käyttää autoonsa enemmän rahaa tulojen noustessa. 



Jos emme ota sukupuolten välistä eroa huomioon, saamme seuraavanlaiset regressioviivat



Vastaava kuva, kun lisäämme interaction-termin malliin:
Interaction-termin lisääminen parantaa mallin ennusteita selkeästi.

lauantai 26. maaliskuuta 2016

R taloustieteellisen tehtävän tukena

Esittelen nyt R:n käyttöä erään taloustieteellisen tehtävän kanssa. R:lla voi tehdä laskutoimituksia ja erityisen hyödyllinen se on kuvien tekemisessä. Tehtävä on seuraavanlainen: 

Hyödykkeen  X

  kysyntä D  on        p = 7 - q/2

  ja tarjonta S on      p = 3/4 + q/2

  määrää hyödykkeen  tasapainohinta ja kysytty  määrä.

Määrää myös suurin mahdollinen veron tuotto, jonka julkinen valta voi 
saada asettamalla hyödykkeelle  t:n suuruisen yksikköveron.

Huom. Käytän * merkkiä tarkoittamaan kertomista, eli 2*2 on kaksi kertaa kaksi. Vastaavasti ^ on potenssimerkki, eli 2^2 tarkoittaa kaksi potenssiin kaksi.

Tässä tehtävässä tulee asettaa p (kysyntä)=p (tarjonta), eli kun piirtää kysyntäkäyrän ja tarjontakäyrän, niin kohdassa p(kysyntä)=p(tarjonta) käyrät leikkaavat (kyseisellä hinnalla kysyntä on sama kuin tarjonta, joten markkinat ovat tasapainossa.

Asetetaan p(kysyntä)=p(tarjonta) ja ratkaistaan määrä q → 7-q/2=3/4+q/2 → q=25/4

Markkinoiden tasapainomäärä on siis 25/4 kappaletta hyödykettä. Markkinahinta saadaan nyt sijoittamalla saatu arvo q=25/4 kysyntäkäyrän- tai tarjontakäyrän yhtälöön (kumpikin tuottaa saman hinnan, koska ratkaisimme määrän q asettamalla ensin hinnat samaksi):

Sijoitetaan tasapainomäärä kysyntäkäyrän yhtälöön p=7-(25/4)*1/2 → p=(56/8)-(25/8) → p=31/8
Vastaus on siis: Hyödykkeen X tasapainohinta: p=31/8 ja Kysytty määrä: q=25/4

R:lla voidaan ratkaista tasapainohinta ja -määrä sekä luoda kuva markkinatasapainosta seuraavasti:
R-Koodi:

kysyntä <- function(q) { 7-(1/2)*q } #Kysyntäkäyrän funktio
plot(kysyntä,from=1,to=10,xlab="Määrä q",ylab="Hinta p",lwd=2) # Kuva kysyntäkäyrästä
tarjonta <- function(q) {3/4+(1/2)*q } # Tarjontakäyrän funktio
curve(tarjonta,add=TRUE,col="blue",xlab="Määrä",lwd=2) # Lisätään kuvaan tarjontakäyrän
legend("topright",c("Kysyntäkäyrä","Tarjontakäyrä"),lty=c(1,1),col=c("black","blue"),lwd=c(2,2))
g <- function(q) { 7-(1/2)*q-3/4-(1/2)*q } # Kysyntä=Tarjonta
uniroot(g,lower=1,upper=10) # Ratkaistaan edellinen funktio.
kysyntä(6.25) # Ratkaistaan tasapainohinta sijoittamalla tasapainomäärä kysyntäkäyrän funktioon
text(7.5,3.9,"(6.25,3.875)") # Lisätään kuvaan tasapainohinnan ja -määrän koordinaatit



Seuraavaksi ratkaistaan suurin mahdollinen veron tuotto, jonka julkinen valta voi saada asettamalla hyödykkeelle t:n suuruisen yksikköveron. Yksikkövero tarkoittaa, että jokaisesta myydystä tuotteesta maksetaan vero. Julkinen valta tietää kysyntäkäyrän (p=7-q/2) sekä tarjontakäyrän (p=3/4+q/2) ja on päättänyt asettaa veron yrityksille. Yritykset laittavat yksikköveron T tarjontakäyrän yhtälöön, eli yhtälö on nyt muotoa: p=3/4+q/2+T
Veron seurauksena tarjontakäyrä nousee siis veron verran ylös, jolloin kysyntä vähenee. Tästä aiheutuu se, ettei markkinahinta nouse veron verran. 

Ratkaistaan ensin tasapainomäärä, kun p(kysyntä)=p(tarjonta+vero). Tämän jälkeen muodostamme verotulojen funktion, jonka derivoimme veron T suhteen. Derivaatta asetetaan nollaksi, jolloin voimme ratkaista verotulot maksimoivan veron määrän. Verotulojen funktio on määrä*vero, eli jos vaihdettu määrä on 8 ja vero 1 niin valtio saa verotuloja 1*8=8. Haluamme siis maksimoida nuo verotulot. Valtio voi nostaa veroa esim. tuosta yhdestä eurosta kahteen euroon, mutta tällöin kysytty määrä vähenee, eli se ei olekaan enää tuo kahdeksan. Se, kuinka paljon kysyntä vähenee, riippuu kysyntäkäyrän muodosta. Kysynnän hintajousto liittyy tähän, eli jos kysyntä laskee hinnan noustessa suhteessa enemmän kuin hinta (joustava kysyntä), niin käytetty rahamäärä vähenee, ja täten myös verotulot tippuvat. Toistaalta, jos kysyntä vähenee hinnan nousteessa suhteessa vähemmän kuin hinta, niin käytetty rahamäärä kasvaa ja verotulot kasvavat (joustamaton kysyntä).


Eli nyt ratkaistaan tehtävä asettamalla p(kysyntä)=p(tarjonta+vero) ja ratkaistaan q →
7-q/2=3/4+q/2+T → q=25/4-T
Verotulojen funktio on qT eli kerrotaan yhtälön kummatkin puolet T:lla ja saadaan → qT=(25/4)T-T^2

Nyt meillä on verotulojen funktio, mutta emme tiedä, mikä veron määrä T maksimoi tämän funktion. Meidän täytyy derivoida verotulojen funktio veron T suhteen. Derivointi siis kertoo meille funktion käyttäytymisestä, eli jos derivaatta on positiivinen, niin itse funktio kasvaa. Derivaatan nollakohdassa derivaatan merkki vaihtuu, joten nämä ovat funktion maksimeja tai minimejä. Jos derivoimme verotulojen funktion T:n suhteen ja emme aseta derivaattaa nollaksi, niin emme tiedä, missä kohtaa funktio saa “ääriarvonsa” eli verotulot maksimoituvat.

Derivoidaan siis qt ja asetetaan derivaatta nollaksi. d(qT)=(25/4)-2T=0
Nyt voimme ratkaista yhtälöstä T:n → 2T=25/4 → T=25/8
Saamamme T on siis maksimiarvo, koska verotulojen funktion toinen derivaatta T:n suhteen on negatiivinen (-2)
Voimme katsoa verotulojen funktiosta verotulojen määrän tällä T:n arvolla. Se on 9,77
(qt=(25/4)*(25/8)-(25/8)^2=9,77)

Voimme myös ratkaista uuden markkinatasapainon hinnan ja määrän. Ratkaistaan ensin määrä q sijoittamalla T=25/8 yhtälöön q=25/4-T
eli yhtälö jonka saimme, kun asetimme p=p. Ratkaistaan: q=25/4-25/8 → q=25/8
Ja nyt voimme tarkistaa tasapainohinnan kysyntä- tai tarjontakäyrän yhtälöstä. Esim. kysyntäkäyrästä: p= 7-q/2 → p=7-(25/8)(1/2) → p=7-25/16 =5,44

Näemme siis, että veron takia uusi tasapainomäärä (25/4 > 25/8) on pienempi ja hinta suurempi 31/8<5,44) kuin ennen veroa.

On myös hyvä havaita, että tuottajien saama hinta on markkinahintaa pienempi, koska he maksavat veron. Tuottajien saama hinta on siis P-T. Hinnan voimme ratkaista tuottajien tarjontakäyrästä → P=3/4+q/2+T → P-T=3/4+q/2
Tasapainossa määrä oli 25/8, joten tuottajien saama hinta on P-T=3/4+25/16=37/16

Alla vielä Kuva uudesta tasapainosta ja kuvan R-koodi.
R-koodi:

kysyntä <- function(q) { 7-(1/2)*q }
plot(kysyntä,from=1,to=10,xlab="Määrä q",ylab="Hinta p",lwd=2)
tarjonta <- function(q) {3/4+(1/2)*q }
curve(tarjonta,add=TRUE,col="blue",xlab="Määrä",lwd=2)
legend(6,6.55,c("Kysyntäkäyrä","Tarjontakäyrä","Tarjontakäyrä(T)"),lty=c(1,1,1),col=c("black",
"blue","red"),lwd=c(2,2,2))
g <- function(q) { 7-(1/2)*q-3/4-(1/2)*q }

uniroot(g,lower=1,upper=10)
kysyntä(6.25)
text(7.5,3.9,"(6.25,3.875)")

g <- function(q) { 7-(1/2)*q-3/4-(1/2)*q-25/8 } #Kysyntä=Tarjonta
uniroot(g,lower=1,upper=10) #Uusi markkinatasapainomäärä
kysyntä(3.125) # Uusi markkinatasapainohinta
tarjontaverolla <- function (q) {3/4+(1/2)*q+25/8 }
curve(tarjontaverolla,add=TRUE, col="red",lwd=2) # Lisätään uusi tarjontakäyrä(T)
text(4.45,5.45,"(3.125,5.438)")