Wstęp do pakietu R

Pakiet R można pobrać i zainstalować ze strony – http://www.r-project.org

Pomoc:

- Jeśli znamy nazwę funkcji, np.

?solve

lub

help(solve)

- Jeśli nie znamy nazwy funkcji:

help.start()

lub

help.search("szukany tekst")

Przypisywanie wartości zmiennym

Zmiennym można przypisywać wartości poprzez operatory przypisania

->

,

<-

,

=

lub poprzez funkcje

assign()

. Na przykład:

a<-5; 4->b; c=b; assign("x",10)
a; b; c; x

Podstawowe struktury danych

-

vector

Do tworzenia wektora zawierającego konkretne elementy używa się funkcji

c()

Na przykład:

a<-c(10,-15,4,-2,8,7)
a

Operator

:

służy do generowania ciągu kolejnych liczb całkowitych, a

funkcja

seq()

do generowania liczb z zakresu z określonym krokiem.

Przykłady:

-10:5
seq(-5,5,by=.2)

Ciąg o powtarzających się elementach tworzy się za pomocą funkcji

rep()

,

np.:

rep(c(1,2),4)

Na wektorach można wykonywać działania lub podstawowe funkcje, np.:

x<-c(2,4,5); y<-c(1,3,2)
x+y; x*y+2
log(x)

Indeksowanie wektorów:

a[2] # drugi element
a[-2] # wszystkie oprócz drugiego
a[c(1,5)] # pierwszy i piaty
a>0 # wektor logiczny - które większe od zera
a[a>0] # elementy większe od zera

Często łatwiej zapamiętać nazwy niż wartości liczbowe. Tekstowe nazwy
przypadków można nadać w np. w następujący sposób:

owoce<-c(1,5,8,10)
names(owoce)<-c("jablko","gruszka","banan","sliwka")
obiad<-owoce[c("jablko","banan")]

-

factor

Faktor jest strukturą przechowującą informacje o powtórzeniach takich
samych wartości, np,:

vect<-c(1,2,1,3,1,2,4,3,2); vect
fact<-factor(vect); fact
summary(fact)
plot(fact)

-

array

Tablice można utworzyć z wektora poprzez określenie wymiarów funkcja

dim()

,

np.:

tbl=1:20
dim(tbl)=c(4,5) # wektor staje się tablica o wymiarach 4,5
tbl
tbl[2,3] # element na pozycji 2,3
tbl[2,] # drugi wiersz
tbl[,3] # trzecia kolumna

Można również dodawać kolumny:

tbl<-cbind(tbl,3); tbl

lub wiersze:

tbl<-rbind(log(tbl[1,]),tbl); tbl

-

data frame

Data frame jest typem tablicy, w której kolumny odpowiadają zmiennym, a
wiersze obserwacjom. Data frame można utworzyć na kilka sposobów:

– wczytując z pliku tekstowego –

read.table("nazwa pliku")

,

– zamieniając tablicę poprzez funkcje

as.data.frame()

,

– poprzez funkcje

data.frame()

, np.:

dat<-data.frame("x"=c(2,4,1,2),"y"=c(1,1,1,1))

– interaktywnie zmieniając istniejącą data.frame za pomocą funkcji

fix()

lub

edit()

.

Do poszczególnych zmiennych odwołujemy się przez

nazwa.obiektu$nazwa.zmiennej

.

Aby do zmiennych można było odwołać się bezpośrednio można dany obiekt
przypisać przez funkcje

attach()

.

Wczytywanie danych z pakietu i z pliku, zapis danych
do pliku tekstowego.

•

Zainstaluj bibliotekę

faraway

, otwórz ja i otwórz dane

pima

z tej

biblioteki:

Z górnego paska wybierz Packages i Instal packages, następnie wybierz
serwer i bibliotekę. Następnie wpisz:

library(faraway); data(pima)

Wszystkie dostępne bazy danych uzyskujemy przez komendę

data()

, a opis

danych poprzez

?pima

. Dane te dotyczą badan nad cukrzyca u Indianek Pima.

Zawierają one ilość ciąży (pregnant), koncentracja glukozy (glucose),
ciśnienie rozkurczowe (diastolic), grubość fałdy skórnej (triceps),
stężenie insuliny (insulin), indeks masy do wzrostu (bmi), funkcja cukrzycy
(diabetes), wiek (ages), czy pacjent wykazuje symptomy cukrzycy (test).

Przeglądnij dane i wyświetl podstawowe statystyki dotyczące tych danych:

summary(pima)

Czy dostrzegasz cos niepokojącego?

W R zmienne brakujace oznacza się przez

NA

. Nanieś niezbędne poprawki:

pima$diastolic[pima$diastolic == 0] <- NA
pima$glucose[pima$glucose == 0] <- NA
pima$triceps[pima$triceps == 0] <- NA
pima$insulin[pima$insulin == 0] <- NA
pima$bmi[pima$bmi == 0] <- NA

Czy zmienna test jest zmienna ilościową czy jakościową? Wskazówka:

pima$test <- factor(pima$test)
levels(pima$test) <- c("negative","positive")

Przeglądnij ponownie podsumowanie danych.

•

Stwórz obiekt

westwood

typu data frame i wypełnij go danymi:

LotSize

ManHours

30

73

20

50

60

128

80

170

40

87

50

108

60

135

30

69

70

148

60

132

Zapisz do pliku ”westwood.txt”:

westwood<-data.frame()
fix(westwood)

Ustaw katalog roboczy z górnego paska.

write.table(westwood,"westwood.txt")

•

Otwórz dane zawierające średnie roczne temperatury z pliku

”temperature.txt”.

Przeglądnij je i zastanów się czy wszystko jest w porządku. Spróbuj wczytać
je poprawnie.

ROZKŁADY ZMIENNYCH

Pakiet R posiada wbudowane algorytmy pozwalające na obliczanie gęstości,
dystrybuanty i kwantyli najczęściej stosowanych rozkładów. Może również
pracować jako precyzyjny generator liczb losowych. Standardowo dostępne są
następujące rozkłady:

beta, binom, cauchy, chisq, exp, f, gamma, geom, hyper, lnorm, logis,
nbinom, norm, pois, t, unif, weibull, wilcox.

Poprzedzając nazwę rozkładu litera d uzyskujemy funkcje gęstości rozkładu.
Analogicznie poprzedzając nazwę literą p uzyskujemy wartości dystrybuanty.
Funkcja kwantylowa (poprzedzona q) podaje taki kwantyl, który po lewej
stronie wykresu gęstości zawiera określone prawdopodobieństwo. Generator
liczb losowych dostępny jest przy poprzedzeniu nazwy litera r.

d

– gęstość rozkładu,

p

– wartość dystrybuanty,

q

– kwantyl,

r

– generator liczb losowych.

Funkcje te pozwalają na traktowanie pakietu R jako zestawu bardzo
dokładnych tablic statystycznych.
Popatrzmy na przykłady:

dnorm(0)

# gęstość rozkładu normalnego w zerze

pnorm(1)-pnorm(-1)

# ile wartości mieści się w N(0,1) w przedziale (-1,1) ?

qt(0.995,5)

# wartość krytyczna t-Studenta dla 99% i 5 stopni swobody

qchisq(0.01,5)

# wartość krytyczna chi-kwadrat dla 5 st. swobody i 99%

dpois(0,0.1)

# wartość prawdopodobieństwa. Poissona dla

1

.

0

=

λ

i n=0

qf(0.99,5,6)

# wartość krytyczna testu F dla n

1

=5, n

2

=6

Kilku słów wymaga wartość 0.995 (nie 0.99) w wywołaniu funkcji rozkładu
t-Studenta. Rozkład ten jest zwykle stosowany w kontekście dwustronnym,
dlatego obszar krytyczny dzielimy równomiernie na obu „końcach” rozkładu.
99% ufności oznacza, ze krytyczny 1% jest podzielony na 2 końce i zawiera
się w przedziałach (0, 0.05) oraz (0.995, 1). Wartość tablicowa jest
kwantylem obliczonym dla takiego właśnie prawdopodobieństwa. Analogicznie
np. dla 95% będzie to 0.975, a dla 99.9% — 0.9995.


Korzystając z funkcji generatora losowego można generować dowolne ciągi
danych do późniejszej analizy. Wygenerujmy przykładowy zestaw 30 liczb o
średniej 50 i odchyleniu standardowym 5, posiadających rozkład normalny:

rnorm(30,50,5)

Oczywiście funkcja ta wygeneruje za każdym razem całkowicie inne wartości,
dlatego tez prowadząc analizy należy je zapamiętać w zmiennej, a potem
używać tej zmiennej w dalszych operacjach.

Warto przy okazji wspomnieć o funkcji sample, generującej wektor danych
wylosowanych z innego wektora.
Np. funkcja

sample(1:6,10,replace=T)

symuluje 10 rzutów kostka (losowanie ze zbioru 1:6), a dane mogą się
powtarzać. Jeśli nie jest podana liczba losowanych danych (np.

sample(1:6)

), funkcja generuje losowa permutacje wektora podanego jako parametr.

WYKRESY

Wykresy słupkowe

Wykresy słupkowe osiągalne są w R poleceniem

barplot

. Popatrzmy na

przykładowe wykresy danych, dotyczących spędzanych godzin przed telewizorem
i komputerem u 100 osób. Wygenerujemy te dane używając funkcji

rnorm

,

zaokrąglając do liczb całkowitych funkcja

round

. Funkcja

abs

zwraca

wartości bezwzględne z danych, gdyż mogą wśród nich trafić się liczby
ujemne.

tv = abs(round(rnorm(100,4,1.5)))
komp = abs(round(rnorm(100,6,2)))
tv
komp

Następnie będziemy produkować różne przykładowe wykresy komendą barplot:

1.

barplot(tv)

—

najprostszy

wykres,

zawierający

100

słupków

odpowiadających poszczególnym danym.
2.

barplot(sort(tv))

— wykres posortowanych danych, pozwalający na wizualna

ocenę ilości poszczególnych odpowiedzi.
3.

barplot(table(tv,tv))

— ilość osób oglądających telewizje godzinę, dwie

etc.
4.

barplot(table(komp,tv),legend.text=T)

— wykres analogiczny, lecz słupki

pokazują

dodatkowo

poszczególne

odpowiedzi

dotyczące

komputera

(są

podzielone). Opcjonalny parametr legend.text pozwala na zamieszczenie
legendy.
5.

barplot(table(komp,komp))

— taki sam wykres, lecz dla czasu spędzonego

przed komputerem.
6.

barplot(table(tv,komp),legend.text=T)

— analogiczny, z podzielonymi

słupkami proporcjonalnie do czasu oglądania telewizji.
7.

barplot(table(komp,tv),beside=T,legend.text=T)

— te same dane, jednak

słupki położone są obok siebie, a nie na sobie.
8.

barplot(table(komp,tv),col=rainbow(10),legend.text=T)

— tutaj listę

tęczowych kolorów generujemy funkcja

rainbow

.

Istnieje kilka funkcji do generowania zestawów kolorów, analogicznie do

rainbow

, np.

heat.colors, terrain.colors, topo.colors, cm.colors

. Istnieje

też możliwość wyspecyfikowania konkretnej listy kolorów, np.

col=c("red","orange","yellow","green","blue",violet")

Wykresy kołowe

Wykresy kołowe generujemy komenda

pie

, z analogiczna składnia. Przykładowo,

poprzez wywołanie

pie(table(tv))

Histogramy

W

czasie

analizy

statystycznej

niejednokrotnie

zachodzi

potrzeba

przedstawienia graficznego rozkładu danych statystycznych. Pakiet R pozwala
na zobrazowanie danych na wykresach w wieloraki sposób.

W pierwszej kolejności załadujemy komendą

data(rivers)

dane dotyczące

długości rzek w USA. Najprostszym sposobem zobrazowania rozkładu tych
danych jest wykres tekstowy stem-and-leaf :

stem(rivers)

Wykres ten nie jest często stosowany, aczkolwiek warto go znać — bardzo
łatwo wkleić go wszędzie tam, gdzie operujemy tylko tekstem. Typowo
„graficzne” wykresy to najczęściej histogramy i tzw. wykresy gęstości
rozkładu. Oto kilka przykładów:

1.

hist(rivers)

— zwykły histogram.

2.

hist(rivers,40)

— tak samo, lecz narzucamy liczbę „przedziałów”

histogramu.
3.

hist(rivers,breaks=c(10,20,50,100,200,300,400,500,1000,4000))

— tutaj

narzucamy konkretne miejsca, w których kończą się przedziały.
4.

plot(density(rivers))

— wykres gęstości.

Niekiedy rozkład danych doświadczalnych prezentuje się na tzw. wykresie
kwantylowo-normalnym (quantile-normal plot). Jest to zależność pomiędzy
wartościami zmiennej, a kwantylami rozkładu normalnego. W idealnym
przypadku, jeśli rozkład zmiennej jest czysto normalny, wykres ten
przedstawia linie prosta. Popatrzmy na trzy przykłady:

1.

qqnorm(rnorm(500))

— wykres dla 500 liczb o rozkładzie normalnym.

2.

qq=rt(500,1);qqnorm(qq);qqline(qq)

— wykres dla 500 liczb o rozkładzie

t-Studenta o 5 stopniach swobody. Funkcja qqline dodaje do wykresu linie
prosta „uśredniającą” rozkład,
3.

qqnorm(runif(500,0,1))

—

wykres

dla

500

liczb,

równomiernie

rozmieszczonych w przedziale (0, 1).

Wykresy dwuwymiarowe

Funkcja plot jest uniwersalnym narzędziem pakietu R, a wykresy generowane
przy jej użyciu mogą być skrajnie odmienne, w zależności od argumentów tej
funkcji. W najprostszym przypadku funkcja ta generuje dwuwymiarowy wykres
punktowy, przedstawiający zależność jednej zmiennej od drugiej. Załadujmy
komenda

data(trees); attach(trees)

dane dotyczące wysokości, średnicy i

objętości pnia 31 drzew, po czym wykonajmy wykresy:

1.

plot(Girth,Volume,type="p")

— wykres punktowy.

2.

plot(Girth,Volume,type="l")

— wykres liniowy.

3.

plot(Girth,Volume,type="b")

— wykres punktowo-liniowy.

4.

plot(Girth,Volume,type="h")

— wykres „kreskowy”.

5.

plot(Girth,Volume,type="s")

— wykres schodkowy.

6.

plot(Girth,Volume,type="S")

— wykres schodkowy (o odwrotnych schodkach).

Wykresy funkcji

R umożliwia rysowanie wykresów dowolnych funkcji matematycznych, przy czym
mogą

one

stanowić

samodzielne

wykresy,

lub

tez

być

dodawane

do

istniejących. Służy do tego funkcja

curve

. Oto kilka przykładów:

1.

curve(pnorm(x),xlim=c(-4,4));curve(pt(x,1),lty=2,add=T)

— rysuje wykres

dystrybuanty

rozkładu

normalnego

linia

ciągła,

a

następnie

dodaje

dystrybuantę rozkładu t-Studenta o 1 stopniu swobody.
2.

hist(rnorm(500),prob=TRUE);curve(dnorm(x),add=T)

— rysuje histogram 500

losowych liczb o rozkładzie normalnym wraz z „idealna” krzywa rozkładu.
3.

curve(sin(x),xlim=c(-2*pi,2*pi)); curve(cos(x),col="red",add=T,xlim=c(-

pi,pi))

—

rysuje wykres funkcji y = sin x w przedziale (−2Π, 2Π), a następnie dodaje
wykres y = cos x, czerwona linia, w przedziale (−Π, Π).

ĆWICZENIA

Zad1.
Zgadnij, a następnie sprawdź jakie będą wyniki następujących poleceń:

> y = c(2, 3, 5, 7, 11, 13)
> length(y)
> x = 2 * 1:5 - 1
> length(x)
> x + y
> sum(x > 5 | x < 3)
> y[-(3:5)]
> y[x]

Zad2.
Twoje czasy dojazdu na uczelnię przez ostatnie dwa tygodnie (10 dni; w
minutach) to:

17 16 20 24 22 15 21 15 17 22

Jakie były najdłuższy, średni i minimalny czasy dojazdu? Jakie było
odchylenie standardowe czasu dojazdu? Ile razy dojazd zajął Ci mniej/więcej
niż średnia -/+ odchylenie standardowe?
Jakie

były

średnie

czasy

dojazdu

dla

wartości

poniżej/ponad

pierwszym/trzecim kwartylem?

Zad3.
Wczytaj wbudowany zbiór danych mtcars.
Zobacz czego dotyczy i sprawdź:

1. ile wynosi maksymalny przebieg (w milach/galon) i który samochód go

osiągnął?

2. jaki wygląda pierwsza trójka samochodów o największej liczbie koni

mechanicznych?

3. jakie są średnie przyspieszenia i odchylenie standardowe liczby koni

dla:

•

wszystkich samochodów,

•

samochodów z / bez automatycznej skrzyni biegów,

•

mercedesów,

Zad4.
Spośród średniej i mediany - estymatorów środka populacji o zadanym
rozkładzie, chcemy wybrać ten o mniejszej wariancji. W tym celu oblicz
stosunek wariancji średniej do wariancji mediany na podstawie n = 1000 prób
dla prób populacji liczności m = 100. Za rozkład populacji przyjmij:
1. Norm(0,1),
2. t(2) (rozkład t z dwoma stopniami swobody),

Więcej na temat systemu R możesz znaleźć w
”An Introduction to R”
dostępnym na stronie projektu R.