BODDy Workshop Data Preparation

library(rgdal)

## Loading required package: sp

## rgdal: version: 1.1-10, (SVN revision 622)
##  Geospatial Data Abstraction Library extensions to R successfully loaded
##  Loaded GDAL runtime: GDAL 1.11.3, released 2015/09/16
##  Path to GDAL shared files: /usr/local/Cellar/gdal/1.11.3_1/share/gdal
##  Loaded PROJ.4 runtime: Rel. 4.9.2, 08 September 2015, [PJ_VERSION: 492]
##  Path to PROJ.4 shared files: (autodetected)
##  Linking to sp version: 1.2-3

library(ggplot2)
library(ggmap)
library(dplyr)

## 
## Attaching package: 'dplyr'

## The following objects are masked from 'package:stats':
## 
##     filter, lag

## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union

library(readr)
library(readxl)
library(lubridate)

## 
## Attaching package: 'lubridate'

## The following object is masked from 'package:base':
## 
##     date

library(zoo)

## 
## Attaching package: 'zoo'

## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric

library(knitr)
library(purrr)

## 
## Attaching package: 'purrr'

## The following object is masked from 'package:dplyr':
## 
##     order_by

Bäume - Excel-Download vom FIS-Broker

baeume = read_excel('strassenbaeume-friedrichshain-kreuzberg.xls') %>% mutate(
  `Baumart deutsch` = factor(`Baumart deutsch`),
  `Baumart botanisch` = factor(`Baumart botanisch`),
  `Gattung` = factor(`Gattung`),
  `Pflanzjahr` = as.numeric(`Pflanzjahr`),
  `Standalter` = as.numeric(`Standalter`),
  `Krone Durchmesser in m` = as.numeric(`Krone Durchmesser in m`),
  `Stamm Umfang in cm` = as.numeric(`Stamm Umfang in cm`),
  `Baumhoehe in m` = as.numeric(`Baumhoehe in m`)
) %>% na.omit() %>%
  group_by(Gattung) %>% mutate(n=n()) %>% ungroup() %>% mutate(r = dense_rank(desc(n))) %>% filter(r<=20) %>% select(-r, -n) %>%
  sample_n(10000) # Keep 10k

## Warning in eval(substitute(expr), envir, enclos): NAs durch Umwandlung
## erzeugt

## Warning in eval(substitute(expr), envir, enclos): NAs durch Umwandlung
## erzeugt

## Warning in eval(substitute(expr), envir, enclos): NAs durch Umwandlung
## erzeugt

## Warning in eval(substitute(expr), envir, enclos): NAs durch Umwandlung
## erzeugt

## Warning in eval(substitute(expr), envir, enclos): NAs durch Umwandlung
## erzeugt

baeume %>% select(`Baumart deutsch`, `Baumart botanisch`, `Gattung`, `Pflanzjahr`, `Standalter`, `Krone Durchmesser in m`, `Stamm Umfang in cm`, `Baumhoehe in m`) %>% summary

##                      Baumart deutsch                  Baumart botanisch
##  LINDE                       :1752   TILIA SPEC.               :1752   
##  KAISER-LINDE                :1361   TILIA INTERMEDIA 'PALLIDA':1361   
##  AHORNBLÄTTRIGE PLATANE      : 658   PLATANUS ACERIFOLIA       : 658   
##  SPITZ-AHORN                 : 581   ACER PLATANOIDES          : 581   
##  WINTER-LINDE                : 577   TILIA CORDATA             : 577   
##  BAUM-HASEL - TÜRKISCHE HASEL: 554   CORYLUS COLURNA           : 554   
##  (Other)                     :4517   (Other)                   :4517   
##       Gattung       Pflanzjahr     Standalter     Krone Durchmesser in m
##  TILIA    :4813   Min.   :1760   Min.   :  1.00   Min.   : 0.010        
##  ACER     :1261   1st Qu.:1966   1st Qu.: 18.00   1st Qu.: 4.000        
##  PLATANUS : 681   Median :1983   Median : 33.00   Median : 6.000        
##  CORYLUS  : 554   Mean   :1980   Mean   : 36.25   Mean   : 6.083        
##  CRATAEGUS: 461   3rd Qu.:1998   3rd Qu.: 50.00   3rd Qu.: 8.000        
##  SORBUS   : 295   Max.   :2015   Max.   :256.00   Max.   :80.000        
##  (Other)  :1935                                                         
##  Stamm Umfang in cm Baumhoehe in m   
##  Min.   :   4.00    Min.   :  3.000  
##  1st Qu.:  52.00    1st Qu.:  7.000  
##  Median :  79.00    Median :  7.000  
##  Mean   :  91.05    Mean   :  9.396  
##  3rd Qu.: 117.00    3rd Qu.: 12.000  
##  Max.   :4773.00    Max.   :400.000  
##

nlevels(baeume$`Baumart botanisch`)

## [1] 141

length(unique(baeume$`Baumart botanisch`))

## [1] 102

nlevels(baeume$`Gattung`)

## [1] 38

length(unique(baeume$`Gattung`))

## [1] 20

Höhe in Meter

ggplot(baeume) + geom_histogram(aes(`Baumhoehe in m`))

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

baeume %>% filter(`Baumhoehe in m` > 30) %>% kable

ID	Baum Nr.	Objektnr.	Objektname	Baumart deutsch	Baumart botanisch	Gattung	ALK Nr.	Straße	Hausnr.	Hausnr. Zusatz	Pflanzjahr	Standalter	Krone Durchmesser in m	Stamm Umfang in cm	Baumhoehe in m	Bezirk
00008100:0016b217	58	GÜRTELST	Gürtelstraße	KAISER-LINDE	TILIA INTERMEDIA ‘PALLIDA’	TILIA	-	Gürtelstraße	-	-	2012	4	1.5	28	400	Friedrichshain-Kreuzberg

ggplot(baeume) + geom_histogram(aes(`Krone Durchmesser in m`))

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

baeume %>% filter(`Krone Durchmesser in m` > 40) %>% kable

ID	Baum Nr.	Objektnr.	Objektname	Baumart deutsch	Baumart botanisch	Gattung	ALK Nr.	Straße	Hausnr.	Hausnr. Zusatz	Pflanzjahr	Standalter	Krone Durchmesser in m	Stamm Umfang in cm	Baumhoehe in m	Bezirk
00008100:0016b6da	36	03608	Ratiborstraße	LINDE	TILIA SPEC.	TILIA	-	Ratiborstraße	-	-	2013	3	80	9	3	Friedrichshain-Kreuzberg

ggplot(baeume) + geom_histogram(aes(`Stamm Umfang in cm`))

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

baeume %>% filter(`Stamm Umfang in cm` > 1000) %>% kable

ID	Baum Nr.	Objektnr.	Objektname	Baumart deutsch	Baumart botanisch	Gattung	ALK Nr.	Straße	Hausnr.	Hausnr. Zusatz	Pflanzjahr	Standalter	Krone Durchmesser in m	Stamm Umfang in cm	Baumhoehe in m	Bezirk
00008100:00168bfc	32/3	LENBACHS	Lenbachstraße	KAISER-LINDE	TILIA INTERMEDIA ‘PALLIDA’	TILIA	-	Lenbachstraße	-	-	2001	15	5.5	4773	7	Friedrichshain-Kreuzberg

13 Meter Stammumfang ist das größte, was ich bei einer Kaiserlinde finden konnte.

Werfe unrealistische und NA raus

# wirft auch NA weg
baeume = baeume %>% filter(
  `Stamm Umfang in cm` < 1000 &
  `Krone Durchmesser in m` < 40 &
  `Baumhoehe in m` < 30)

Blick via scatterplots

ggplot(baeume %>% filter(`Krone Durchmesser in m` < 40 & `Baumhoehe in m` < 30)) + geom_point(aes(Standalter, `Krone Durchmesser in m`, color=`Gattung`)) + scale_colour_discrete(guide=F)

baeume %>% filter(Standalter > 200) %>% kable

ID	Baum Nr.	Objektnr.	Objektname	Baumart deutsch	Baumart botanisch	Gattung	ALK Nr.	Straße	Hausnr.	Hausnr. Zusatz	Pflanzjahr	Standalter	Krone Durchmesser in m	Stamm Umfang in cm	Baumhoehe in m	Bezirk
00008100:00163dcf	24	04737	Waterloo-Ufer	SCHEIN-AKAZIE - SILBERREGEN	ROBINIA PSEUDOACACIA	ROBINIA	-	Waterloo-Ufer	-	-	1760	256	5	153	7	Friedrichshain-Kreuzberg

ggplot(baeume %>% filter(`Krone Durchmesser in m` < 40 & `Baumhoehe in m` < 30)) + geom_point(aes(Standalter, `Baumhoehe in m`, color=`Gattung`, size=`Stamm Umfang in cm`)) + scale_colour_discrete(guide=F)

Wind/Wetter

Die Windwerte holen wir uns von DWD CDC (http://www.dwd.de/DE/klimaumwelt/cdc/cdc_node.html). Normalerweise würden wir alle Berliner Wetterstationen haben wollen, aber der einfachheit halber nehmen wir nur die Daten der Station Schönefeld.

wetter = read_delim('produkt_klima_Tageswerte_20141202_20160603_00427.txt',delim=';')

wetter %>% summary

##   STATIONS_ID    MESS_DATUM       QUALITAETS_NIVEAU LUFTTEMPERATUR   
##  Min.   :427   Min.   :20141202   Min.   :1.000     Min.   :-10.200  
##  1st Qu.:427   1st Qu.:20150418   1st Qu.:3.000     1st Qu.:  3.700  
##  Median :427   Median :20150902   Median :3.000     Median :  8.350  
##  Mean   :427   Mean   :20152873   Mean   :2.891     Mean   :  9.168  
##  3rd Qu.:427   3rd Qu.:20160118   3rd Qu.:3.000     3rd Qu.: 14.175  
##  Max.   :427   Max.   :20160603   Max.   :3.000     Max.   : 29.400  
##    DAMPFDRUCK     BEDECKUNGSGRAD  LUFTDRUCK_STATIONSHOEHE  REL_FEUCHTE   
##  Min.   : 2.200   Min.   :0.000   Min.   : 968.3          Min.   :38.58  
##  1st Qu.: 6.500   1st Qu.:4.300   1st Qu.:1004.6          1st Qu.:65.73  
##  Median : 8.200   Median :5.700   Median :1010.6          Median :78.81  
##  Mean   : 9.121   Mean   :5.411   Mean   :1010.0          Mean   :76.45  
##  3rd Qu.:10.900   3rd Qu.:6.900   3rd Qu.:1016.0          3rd Qu.:87.83  
##  Max.   :21.700   Max.   :8.000   Max.   :1032.0          Max.   :99.92  
##  WINDGESCHWINDIGKEIT LUFTTEMPERATUR_MAXIMUM LUFTTEMPERATUR_MINIMUM
##  Min.   :-999.000    Min.   :-8.200         Min.   :-12.900       
##  1st Qu.:   3.000    1st Qu.: 6.925         1st Qu.:  0.125       
##  Median :   4.100    Median :12.300         Median :  3.900       
##  Mean   :  -4.665    Mean   :13.499         Mean   :  4.651       
##  3rd Qu.:   5.300    3rd Qu.:19.900         3rd Qu.:  8.975       
##  Max.   :  11.400    Max.   :37.800         Max.   : 21.000       
##  LUFTTEMP_AM_ERDB_MINIMUM WINDSPITZE_MAXIMUM  NIEDERSCHLAGSHOEHE
##  Min.   :-17.900          Min.   :-999.0000   Min.   : 0.000    
##  1st Qu.: -2.400          1st Qu.:   8.2000   1st Qu.: 0.000    
##  Median :  1.000          Median :  10.6000   Median : 0.000    
##  Mean   :  2.012          Mean   :   0.4449   Mean   : 1.245    
##  3rd Qu.:  6.300          3rd Qu.:  13.9000   3rd Qu.: 1.100    
##  Max.   : 20.400          Max.   :  29.6000   Max.   :23.800    
##  NIEDERSCHLAGSHOEHE_IND SONNENSCHEINDAUER  SCHNEEHOEHE     
##  Min.   :0.000          Min.   : 0.000    Min.   : 0.0000  
##  1st Qu.:0.000          1st Qu.: 0.283    1st Qu.: 0.0000  
##  Median :6.000          Median : 3.733    Median : 0.0000  
##  Mean   :4.253          Mean   : 4.643    Mean   : 0.1855  
##  3rd Qu.:6.000          3rd Qu.: 7.929    3rd Qu.: 0.0000  
##  Max.   :8.000          Max.   :15.683    Max.   :10.0000  
##      eor           
##  Length:550        
##  Class :character  
##  Mode  :character  
##                    
##                    
##

Fehlende Werte sind als -999 kodiert. Das kommt sowohl beim Wind als auch bei den Windspitzen vor, die uns interessieren.

Wir ersetzen die -999er Werte durch NA und dann füllen wir mit dem letzen vorherigen Wert auf.

wetter = wetter %>%
  transmute(Datum=ymd(MESS_DATUM), Windmittel=WINDGESCHWINDIGKEIT, Windspitze=WINDSPITZE_MAXIMUM) %>%
  mutate(Windmittel=na.locf(replace(Windmittel, Windmittel==-999, NA)), Windspitze=na.locf(replace(Windspitze, Windspitze==-999, NA)))

ggplot(wetter) + geom_line(aes(Datum, Windmittel)) + geom_line(aes(Datum, Windspitze), color='red')

ggplot(wetter) + geom_histogram(aes(Windspitze))

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

Ab 20 m/s gilt als Sturm.

max(wetter$Datum)

## [1] "2016-06-03"

Zufallsinspektionen (Unbenutzt!)

Nun generieren wir uns Zufallsinspektionen. Laut einem Ordnungsamt-Online-Eintrag gilt:

Der öffentliche Baumbestand wird ein mal jährlich hinsichtlich der Gewährleistung der Verkehrssicherkeit kontrolliert.

inspections = baeume %>% transmute(
  ID=ID,
  Datum=sample(seq(max(wetter$Datum) - years(1), max(wetter$Datum), by=1), n(), replace=T),
  Risiko=sample(seq(0,3), n(), replace=T, prob=c(.9,.09,.009,.001))
  )

Zufalls-Aufräumarbeiten

Gedanken

Welche Baumarten sind bei Sturm gefährdet? Nach internetrecherche sind

Flachwurzler

PSEUDOTSUGA (Fichte)
SALIX (Weide)

Herzwurzler

ACER (Ahorn)
CARPINUS (Buche)
LARIX (Lärche)
TILIA (Linde)
CORYLUS (Baumhasel)
LIQUIDAMBAR (Amberbaum)
PRUNUS (Kirsche)
SORBUS (Mehlbeere)
FRAXINUS (Esche)
BETULA (Birke)

Pfahlwurzler

PINUS (Kiefer)
QUERCUS (Eiche)
ROBINIA (Scheinakazie) … alle

FLACHWURZLER_RISIKO = 1
HERZWURZLER_RISIKO = 0.5
PFAHLWURZLER_RISIKO = 0.25
risiko_gattung = function(x) {
  x = c(
    PSEUDOTSUGA=FLACHWURZLER_RISIKO,
    SALIX=FLACHWURZLER_RISIKO,
    ACER=HERZWURZLER_RISIKO,
    CARPINUS=HERZWURZLER_RISIKO,
    LARIX=HERZWURZLER_RISIKO,
    TILIA=HERZWURZLER_RISIKO,
    CORYLUS=HERZWURZLER_RISIKO,
    LIQUIDAMBAR=HERZWURZLER_RISIKO,
    PRUNUS=HERZWURZLER_RISIKO,
    SORBUS=HERZWURZLER_RISIKO,
    FRAXINUS=HERZWURZLER_RISIKO,
    BETULA=HERZWURZLER_RISIKO
  )[x]
  x[is.na(x)] = PFAHLWURZLER_RISIKO
  return(x)
}

risiko_gattung(c('PSEUDOTSUGA', 'PSEUDOTSUGA', 'BETULA', 'QUERCUS'))

## PSEUDOTSUGA PSEUDOTSUGA      BETULA        <NA> 
##        1.00        1.00        0.50        0.25

# baeume %>% group_by(`Baumart deutsch`) %>% summarise(count=n(), Gattung=first(Gattung)) %>% arrange(count) %>% head()

Zudem sind Risikofaktoren:

Höhe ~ Risiko (ab 20m stark erhöht) Stammdurchmesser ~ 1/Risiko Kronendurchmesser ~ Risiko Risikofaktor ~ Risiko Alter ~ Risiko Wurzeltyp ~ Risiko

Windspitze ~ Risiko Tage seit Schaden ~ Risiko

=> (HöheKronendurchmesser/Stammdicke)Risikofaktor^2

Andere Features (die ich nicht habe):

Wie offen steht der Baum?
Kronenprozent

Daten generieren

max_wind = 25 # 26.5
wetter %>% filter(Windspitze > max_wind)

## Source: local data frame [5 x 3]
## 
##        Datum Windmittel Windspitze
##       (date)      (dbl)      (dbl)
## 1 2015-01-10        9.6       29.2
## 2 2015-03-29        6.7       26.7
## 3 2015-03-30       11.4       26.4
## 4 2015-03-31       11.4       29.6
## 5 2015-11-29        9.1       26.2

Generiere einen Data Frame mit jedem Baum für jedes Sturmereignis

candidates = baeume %>%
  select(ID, Gattung, `Standalter`, `Krone Durchmesser in m`, `Stamm Umfang in cm`, `Baumhoehe in m`) %>%
  mutate(join=1) %>%
  full_join(wetter %>% filter(Windspitze > max_wind) %>% mutate(join=1)) %>% select(-join)

## Joining by: "join"

candidates %>% summary

##       ID                 Gattung        Standalter    
##  Length:49985       TILIA    :24050   Min.   :  1.00  
##  Class :character   ACER     : 6305   1st Qu.: 18.00  
##  Mode  :character   PLATANUS : 3405   Median : 33.00  
##                     CORYLUS  : 2770   Mean   : 36.26  
##                     CRATAEGUS: 2305   3rd Qu.: 50.00  
##                     SORBUS   : 1475   Max.   :256.00  
##                     (Other)  : 9675                   
##  Krone Durchmesser in m Stamm Umfang in cm Baumhoehe in m  
##  Min.   : 0.010         Min.   :  4.0      Min.   : 3.000  
##  1st Qu.: 4.000         1st Qu.: 52.0      1st Qu.: 7.000  
##  Median : 6.000         Median : 79.0      Median : 7.000  
##  Mean   : 6.076         Mean   : 90.6      Mean   : 9.357  
##  3rd Qu.: 8.000         3rd Qu.:117.0      3rd Qu.:12.000  
##  Max.   :30.000         Max.   :565.0      Max.   :25.000  
##                                                            
##      Datum              Windmittel      Windspitze   
##  Min.   :2015-01-10   Min.   : 6.70   Min.   :26.20  
##  1st Qu.:2015-03-29   1st Qu.: 9.10   1st Qu.:26.40  
##  Median :2015-03-30   Median : 9.60   Median :26.70  
##  Mean   :2015-05-02   Mean   : 9.64   Mean   :27.62  
##  3rd Qu.:2015-03-31   3rd Qu.:11.40   3rd Qu.:29.20  
##  Max.   :2015-11-29   Max.   :11.40   Max.   :29.60  
##

stdize = function(x, ...) {(x - min(x, ...)) / (max(x, ...) - min(x, ...))}

candidates = candidates %>% mutate(
  #Höhe ~ Risiko (ab 20m stark erhöht)
  Hoehenrisiko=stdize(`Baumhoehe in m`)*0.99+0.01,
  #Stammdurchmesser ~ 1/Risiko
  Stammrisiko=stdize(-`Stamm Umfang in cm`)*0.99+0.01,
  #Kronendurchmesser ~ Risiko
  Kronenrisiko=stdize(`Krone Durchmesser in m`)*0.99+0.01,
  #Alter ~ Risiko
  Altersrisiko=stdize(Standalter)*0.99+0.01,
  #Wurzeltyp ~ Risiko
  Gattungsrisiko=risiko_gattung(Gattung),
  #Windspitze ~ Risiko
  Windrisiko=stdize(Windspitze)*0.99+0.01,
  #Tage seit Schaden ~ Risiko
  #Risikofaktor ~ Risiko
  Risiko = (Hoehenrisiko*Stammrisiko*Kronenrisiko*Altersrisiko*Gattungsrisiko*Windrisiko)^0.8
 )

ggplot(candidates) + geom_histogram(aes(Risiko))

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

candidates = candidates %>% mutate(Schaden=candidates$Risiko %>% map_int(function(p) sample(seq(0,3), 1, prob=c((1-p)*6/8,(1-p)*1/8,(1-p)/8,p)))) #~ rbinom(1, 3, p=.)))
factor(candidates$Schaden) %>% summary

##     0     1     2     3 
## 37290  6233  6194   268

candidates %>% filter(Schaden > 2) %>% ggplot + geom_point(aes(Risiko, `Standalter`, size=`Stamm Umfang in cm`, color=Gattung))

candidates %>% filter(Schaden > 0) %>% .$Gattung %>% unique() %>% length

## [1] 20

candidates %>% .$Gattung %>% unique() %>% length

## [1] 20

ggplot(candidates) + geom_boxplot(aes(factor(Schaden), Risiko))

Tabelle für Schäden

Schaeden = candidates %>% filter(Schaden > 0) %>% mutate(
  Schadensart = factor(c(NA, 'niedrig (Abgebrochene Zweige)', 'Mittel (Abgebrochene Äste)', 'Hoch (Windwurf/Windbruch)')[Schaden + 1])
) %>% select(ID, Datum, Schadensart)

Schaeden %>% write.csv('schaeden.csv')

Clean Bäume

write.csv(baeume, 'baeume-clean.csv')

Test der logistischen Regression auf den candidates

split_date = ymd('2015-03-31')

training = candidates %>%
  filter(Schaden > 0 & Datum < split_date) %>%
  mutate(high_risk=as.numeric(Schaden==3))

test = candidates %>%
  filter(Schaden > 0 & Datum >= split_date) %>%
  mutate(high_risk=as.numeric(Schaden==3))

model <- glm(high_risk ~ `Gattung` + `Standalter` + `Krone Durchmesser in m` + `Stamm Umfang in cm` + `Baumhoehe in m` + Windspitze, family=binomial(link='logit'),data=training)

p <- predict(model, type="response", newdata = test)

summary(model)

## 
## Call:
## glm(formula = high_risk ~ Gattung + Standalter + `Krone Durchmesser in m` + 
##     `Stamm Umfang in cm` + `Baumhoehe in m` + Windspitze, family = binomial(link = "logit"), 
##     data = training)
## 
## Deviance Residuals: 
##     Min       1Q   Median       3Q      Max  
## -1.2921  -0.1847  -0.1083  -0.0625   3.5251  
## 
## Coefficients:
##                            Estimate Std. Error z value Pr(>|z|)    
## (Intercept)              -2.662e+01  2.261e+00 -11.774  < 2e-16 ***
## GattungAESCULUS          -1.129e+00  6.344e-01  -1.779  0.07523 .  
## GattungBETULA             1.385e-01  1.057e+00   0.131  0.89582    
## GattungCARPINUS          -1.564e+01  1.156e+03  -0.014  0.98920    
## GattungCORYLUS           -5.461e-01  4.588e-01  -1.190  0.23390    
## GattungCRATAEGUS         -1.529e+01  5.398e+02  -0.028  0.97740    
## GattungFRAXINUS          -8.342e-01  7.495e-01  -1.113  0.26571    
## GattungGINKGO            -1.600e+01  1.399e+03  -0.011  0.99088    
## GattungGLEDITSIA         -7.271e-01  6.161e-01  -1.180  0.23792    
## GattungPLATANUS          -9.114e-01  3.083e-01  -2.956  0.00311 ** 
## GattungPOPULUS           -2.178e+00  8.075e-01  -2.697  0.00700 ** 
## GattungPRUNUS            -1.548e+01  8.997e+02  -0.017  0.98627    
## GattungPYRUS             -1.500e+01  7.793e+02  -0.019  0.98464    
## GattungQUERCUS           -1.024e+00  5.577e-01  -1.837  0.06623 .  
## GattungROBINIA           -3.880e-01  7.888e-01  -0.492  0.62277    
## GattungSALIX             -1.655e+01  2.367e+03  -0.007  0.99442    
## GattungSOPHORA           -4.466e-01  1.047e+00  -0.427  0.66962    
## GattungSORBUS            -1.120e-01  6.266e-01  -0.179  0.85818    
## GattungTILIA             -1.068e+00  2.247e-01  -4.754 2.00e-06 ***
## GattungULMUS             -1.737e+01  1.759e+03  -0.010  0.99212    
## Standalter                1.277e-02  6.829e-03   1.870  0.06147 .  
## `Krone Durchmesser in m`  1.038e-01  3.827e-02   2.714  0.00665 ** 
## `Stamm Umfang in cm`     -4.746e-03  3.488e-03  -1.361  0.17355    
## `Baumhoehe in m`          1.317e-01  2.957e-02   4.452 8.49e-06 ***
## Windspitze                7.558e-01  7.770e-02   9.727  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 1416.8  on 7527  degrees of freedom
## Residual deviance: 1129.3  on 7503  degrees of freedom
## AIC: 1179.3
## 
## Number of Fisher Scoring iterations: 18

plot(test$Risiko, p)

library(ROCR)

## Loading required package: gplots

## 
## Attaching package: 'gplots'

## The following object is masked from 'package:stats':
## 
##     lowess

pr <- prediction(p, test$high_risk)
prf <- performance(pr, measure = "tpr", x.measure = "fpr")
plot(prf)

auc <- performance(pr, measure = "auc")
auc <- auc@y.values[[1]]
auc

## [1] 0.8735065

test %>% mutate(prediction = p) %>% arrange(desc(prediction)) %>%
  mutate(high_risk_found=cumsum(high_risk)/sum(high_risk), percentage_checked=row_number()/n()) %>% 
  ggplot() + geom_step(aes(percentage_checked, high_risk_found), color='red') + geom_line(aes(percentage_checked, percentage_checked))