Cviceni VPFMP LS 2012/13

Jan Hurt

20130227 Cviceni 1

FinancialData["USD/CZK"]

Casto pouzivejte Ctrl+S !!!

Vyreste rovnici a (x^2)+b x+c=0.

reseni=Solve[a (x^2)+b x+c==0,x] (* nutne mezery mezi a a x *)

Dosazovani (lokalni pravidlo)

reseni/. {a->5,b->4,c->2}
reseni/. {a->1,b->8,c->4}

a=5
reseni
reseni/. {a->4,b->8,c->2}

Clear[a] (* totez, co a=. *)
reseni

sest=6;

reseni/.{a->sest,b->8,c->2}

x

?Global`*

koreny=x/.reseni

koreny/. {a->4,b->8,c->2}

Pozor!

x=5;
Solve[a (x^2)+b x+c==0,x]

x=.;Solve[a (x^2)+b x+c==0,x]

koreny

1. koren, 2. koren

koreny[[1]]
First[koreny]
koreny[[2]]
Last[koreny]

{koreny[[1]],First[koreny],koreny[[2]],Last[koreny]}

Dvojitou hranatou zavorku [[ vytvorime pomoci Esc[[Esc.

20130306 Cviceni 2

FAQ

Promenne v jedne seanci jsou platne ve vsech soucasne otevrenych noteboocich!!!

abc=5;

Superfunkce FinancialData

FinancialData["^DJI",{2013,2,25}]

Date[]

FinancialData["^DJI",{2013,2,25}]

Pro dany rozsah dat:

FinancialData["^DJI",{{2013,1,25},{2013,2,25}}]

DateListPlot[FinancialData["^DJI",{2003,2,25}]]

{FinancialData["^DJI"],FinancialData["^SP500"]}

Mozna nefunguje pod Windows XP:

membersDJI=FinancialData["^DJI","Members"]

Map[{#,FinancialData[#,"Name"]}&,membersDJI]//TableForm

Definice a grafy funkci

Definujte funkci f(x)=sin(x) na intervalu (0,π], na intervalu (π,10π]  periodicky. Vyuzijte funkci Mod. Graficky znazornete.

Reseni: Funkce Mod je zbytek po deleni:

{Mod[5,3], Mod[(3/2)π,π]}

mujsin[x_]:=Sin[Mod[x,π]]
Plot[mujsin[x],{x,0,10π}]

Alternativne:

Plot[Abs[Sin[x]],{x,0,10π}]

Uloha: V textovem modu jmeno, tucne, zvetsit, na radku vpravo:

Jan Hurt

Symbolicke operace 1

e1=Expand[(((x+y))^7)]
e2=Factor[e1]

Definice polynomu 1+x+(x^2)+...+(x^n)

Vyuzijte funkci Table[x,{n+1}], Range[n], Range[n-1]

Total[{x,y,z}]
Plus[x,y,z]

Clear[poly];
poly[x_,n_]:=Total[((Table[x,{n+1}])^(Range[n+1]-1))]

poly[x,3]

Alternativne

(x^({1,2}))

Clear[poly2];
poly2[x_,n_]:=Total[(x^(Range[n+1]-1))]
poly2[y,3]

e3=poly[e1,3]

e4=e3//Expand

e4//Simplify

poly[0,2]

Dodefinovani (0^0)

Vnitrni zobrazeni mocniny v Mathematica

FullForm[(x^y)]

(0^0)=1

Unprotect[Power];
(0^0)=1
Protect[Power];
(0^0)

poly[0,2]

poly[🐺,3]

20130313 Cviceni 3

FAQ

Soustava linearnich rovnic

Clear[a11,a12,b1,a21,a22,b2];
reseni=Solve[soustava={a11 x1+a12 x2==b1,a21 x1+a22 x2==b2},{x1,x2}]

Reduce[soustava,{x1,x2}]

Symbolicke operace 2

Vytvorujici funkce pravdepodobnosti 1

Jestlize nahodna velicina X nabyva s nenulovou pravdepodobnosti nezapornych celych hodnot s pravdepodobnostmi (Subscript[p, 0]),(Subscript[p, 1]),..., pak G(s)=(∑_(k=0)^∞)(s^k)(Subscript[p, k])=E,(s^X) se nazyva vytvorujici funkce pravdepodobnosti. Zname-li vytvorujici funkci, muzeme pravdepodobnosti (Subscript[p, k]) spocitat pomoci

(Subscript[(Subscript[p, k])=(1/(k!))((d^k)/(d,(s^k)))G(s)|, s=0]).

1. Vyzkousejte na binomickem rozdeleni.

?*Distribution

BinomialDistribution

Clear[n,p];
pdf=PDF[BinomialDistribution[n,p],k]

(* vytvorujici funkce *)
gbin[n_,p_,s_]:=Sum[(s^k)PDF[BinomialDistribution[n,p],k],{k,0,n}]

gbin[n,p,s]

Assuming[n>0,gbin[n,p,s]]

nebo rovnou definovat gbin ... /; n>0.

Jina moznost: G(s)={0,1,...,n}.{(Subscript[p, 0]),...,(Subscript[p, n])}:

gbin2[n_,p_,s_]:=(s^(Range[n+1]-1)).Table[PDF[BinomialDistribution[n,p],k],{k,0,n}]/;n>0
gbin2[10,p,s]
gbin2[10,p,s]//Simplify

Zpetny vypocet pravdepodobnosti (Subscript[p, k]) pro Bi(n,p).

psti[n_,p_]:=Table[(1/(k!))D[gbin[n,p,s],{s,k}],{k,0,n}]/.s->0
psti[10,p]

liplot[n_,p_]:=ListPlot[{Range[n+1]-1,psti[n,p]}//Transpose,PlotRange->All]
GraphicsRow[ {liplot[10,0.5],liplot[100,0.5],liplot[100,0.2]}]

psti[100,0.2][[80]]

20130320 Cviceni 4

FAQ

Vytvorujici funkce pravdepodobnosti 2

2a. (Domaci cviceni) Vytvorte (efektivnim zpusobem) funkci (Subscript[g, 1])(z,p,r):=((p^r)/((1-(1-p)z)(1-(1-p)(z^2))...(1-(1-p)(z^r)))). Vyuzijte vysledku posledniho cviceni 20130306.
2b. (Obtiznejsi, bonifikace +1 stupen u zkousky, kdo mi jako prvni zasle spravne reseni e-mailem). Vyse uvedena funkce (Subscript[g, 1]) neni obecne ve vsech pripadech vytvorujici funkci pravdepodobnosti. Upravte ji tak, aby byla.

Bernoulliovske pokusy 1

Uvod

Posloupnost tzv. bernoulliovskych pokusu je posloupnost nezavislych nahodnych pokusu se dvema moznymi vysledky Zdar resp. Nezdar (zkracene Z resp. N). Vystupem je napr. posloupnost typu ZZNZNNNZZZZZNZZ. Nejjednodussi model je hazeni minci,  Z ztotoznen s licem (drive panna) (L), N ztotoznen s rubem (drive orel, pozdeji lev) (N). Predpoklada se, ze pravdepodobnosti (psti) Z resp. N jsou P(Z):=p) resp. P(N):=1-p=:q. V prikladu hazeni (regulerni, spravedlivou) minci je p=q=1/2. V pripade  hazeni (regulerni, spravedlivou) kostkou muzeme definovat Zdar jako padnuti 6, ostatni pripady definovat jako Nezdar. Zde p=1/6 a q=5/6.

Pro znazsi kvantifikaci se zavadeji nahodne veliciny zvane indikatory. Jakmile narazime na pojem “pocet” (number of, nekdy znacene #No), je dobre si uvedomit toto: koname-li n bernoulliovskych pokusu, pak definujeme nahodne veliciny (Subscript[I, i])=1, jestlize i-ty pokus skoncil zdarem, a (Subscript[I, i])=0, jestlize pokus skoncil nezdarem, i=1,...n. Indikator se take nazyva alternativni nahodna velicina.

Pocet Zdaru (Subscript[X, n])=(∑_(i=1)^n)(Subscript[I, i])

Rozdeleni alternativni nahodne veliciny: BernoulliDistribution

PDF[BernoulliDistribution[p],k]

SeedRandom[27];
RandomChoice[{"Z","N"},20]

Viz tez RandomSample

Clear[nahodneZdaryaNezdary];
nahodneZdaryaNezdary [n_Integer(* pocet pokusu *)]:=RandomChoice[{"Z","N"},n]
nahodneZdaryaNezdary [20]
nahodneZdaryaNezdary [2.1]

Uloha: Pridejte k teto procedure dalsi  volitelny parametr (SeedRandom).

{SeedRandom[23],{5,4}}
{SeedRandom[23],{5,4}}[[2]]

Clear[nahodneZdaryaNezdary2];nahodneZdaryaNezdary2[n_Integer,seed___Integer]:=
{SeedRandom[seed],RandomChoice[{"Z","N"},n]}//Last
 nahodneZdaryaNezdary2[5]
 nahodneZdaryaNezdary2[5,13]
nahodneZdaryaNezdary2[5]
 nahodneZdaryaNezdary2[5,13]

Pozor! Prirazeni musi byt odlozene!!!

nahodneZdaryaNezdarySPATNE [n_Integer(* pocet pokusu *)]=RandomChoice[{"Z","N"},20]
nahodneZdaryaNezdarySPATNE [20]

Jinak by

Module

Clear[nahodneZdaryaNezdary3];nahodneZdaryaNezdary3[n_Integer,seed___Integer]:=
Module[{vnitrnipromennamodulu1,vnitrnipromennamodulu2},
If[Length[{seed}]>0,SeedRandom[seed]];
RandomChoice[{"Z","N"},n]
            ]
 nahodneZdaryaNezdary3[5]
 nahodneZdaryaNezdary3[5,13]
nahodneZdaryaNezdary3[5]
 nahodneZdaryaNezdary3[5,13]
nahodneZdaryaNezdary3[5]
nahodneZdaryaNezdary3[5]

Table[((Count[nahodneZdaryaNezdary3[10000],"Z"])/10000),{10}]//N

Geometricke rozdeleni

Doba cekani na 1. zdar =: K ma geometricke rozdeleni. Jinymi slovy je K pocet nezdaru predchazejici prvnimu zdaru.

pst[p_,k_]:=p(((1-p))^k)
(∑_(k=0)^∞)pst[p,k]

PDF[gp=GeometricDistribution[p],k]

{Mean[gp],Variance[gp]}//TraditionalForm

Vytvorujici funkce (okamzite prirazeni!!!)

Clear[vytvorujicifunkcegeometrickehorozdeleni];vytvorujicifunkcegeometrickehorozdeleni[z_,p_]=(∑_(k=0)^∞)pst[p,k] (z^k)

vytvorujicifunkcegeometrickehorozdeleni[1,p]

20130327 Cviceni 5

FAQ

Bernoulliovske pokusy 2 (Zaporne binomicke rozdeleni)

Necht (Subscript[K, r]) je nahodna velicina s nasledujicim vyznamem: r-ty zdar v posloupnosti bernoulliovskych pokusu nastane v pokusu ((Subscript[K, r])+r)-tem. Je tedy (Subscript[K, r]) pocet nezdaru predchazejicich r-temu zdaru. Rozdeleni nahodne veliciny (Subscript[K, r]) je zaporne binomicke radu r. (nebo tez negativne binomicke). Zrejme (Subscript[K, 1]) ma geometricke rozdeleni a soucet r nezavislych stejne rozdelenych nahodnych velicin s geometrickym rozdelenim ma rozdeleni jako (Subscript[K, r]).

(* Zabudovana funkce *)NegativeBinomialDistribution

Cviceni 1

1.  Spoctete vytvorujici funkci zaporne binomickeho rozdeleni. (Znate minimalne dva zpusoby.)

1. 1 Soucet nezavislych geometrickych rozdeleni

g[s_,p_]=(∑_(k=0)^∞)PDF[GeometricDistribution[p],k] (s^k)

gneg[s_,r_,p_]=((g[s,p])^r)

2. 2 Primo z definice vytvorujici funkce

(∑_(k=0)^∞)PDF[NegativeBinomialDistribution[r,p],k] (s^k)

Cviceni 2

2. Vypoctete odpovidajici pravdepodobnosti zaporne binomickeho rozdeleni, pomoci funkce Timing zjistete casovou narocnost. Zkuste se zamyslet nad tim, proc pro vyssi hodnoty r a k vypocet trva relativne dlouho. Vypoctete stredni hodnotu, rozptyl a smerodatnou odchylku.

pneg[k_,r_,p_]:=(1/(k!))D[gneg[s,r,p],{s,k}]/.s->0

Timing[tneg100=Table[pneg[k,30,1/6],{k,0,1000}];]
N[Total[tneg100],30]

ciselnecharakteristikynegativnebinomickehorozdeleni={Mean[dn=NegativeBinomialDistribution[r,p]],Variance[dn],StandardDeviation[dn]}

ciselnecharakteristikynegativnebinomickehorozdeleni/.{p->1/2,r->10}
ciselnecharakteristikynegativnebinomickehorozdeleni/.{p->1/6,r->10}

Cviceni 3

3. Graficky znazornete pravdepodobnosti pro p=1/6, r=10 (doba cekani na sesty vyskyt sestky minus 5). Pro obecna p a r pouzijte funkci Manipulate.

Zaporne binomicke rozdeleni (reseni)

PDF[NegativeBinomialDistribution[r,p],k]

1.1 Jako soucet geometrickych:

.g[z_,p_,r_]:=(((((p )/(1-(1-p)z))))^r)

(t7=Table[{k,D[((g[z,p,r])/(k!)),{z,k}]/.z->0},{k,0,15}]//Simplify)//TableForm

tot7=Total[Transpose[t7]//Last]//FullSimplify

tot7/.{p->2/3,r->4}

1.2 Primo z definice:

Clear[sng];
sng[z_,p_,r_]=(∑_(k=0)^∞)(z^k)PDF[NegativeBinomialDistribution[r,p],k]

(t8=Table[{k,D[((sng[z,p,r])/(k!)),{z,k}]/.z->0},{k,0,15}]//Simplify)//TableForm

2.

{Mean[NegativeBinomialDistribution[r,p]],StandardDeviation[NegativeBinomialDistribution[r,p]]}

3. Graficke znazorneni zaporne binomickeho rozdeleni

PointSize

Manipulate[ListPlot[{Range[Max[0,Floor[(((1-p) r)/p)-5((Sqrt[(1-p) r])/p)]],Floor[(((1-p) r)/p)+5((Sqrt[(1-p) r])/p)],1],Table[PDF[NegativeBinomialDistribution[r,p],k],{k,Max[0,Floor[(((1-p) r)/p)-5((Sqrt[(1-p) r])/p)]],Floor[(((1-p) r)/p)+5((Sqrt[(1-p) r])/p)]}]}//Transpose,PlotStyle->PointSize[0.012]],{{p,0.5},0.1,0.9},{r,1,10,1}]

Vysvetlivky

Probability[Abs[x-μ]>3σ,x\[Distributed]NormalDistribution[μ,σ]]//N

≈

Soucet (ne nutne stejne rozdelenych) nezavislych n. v.

Jsou-li (Subscript[X, 1]),...,(Subscript[X, n])  iid s vytvorujici funkci G(s), ma nahodna velicina (Subscript[X, 1])+...+(Subscript[X, n]) vytvorujici funkci (G^n)(s). Jsou-li (Subscript[X, 1]),...,(Subscript[X, n]) nezavisle (ne nutne stejne rozdelene!) s vytvorujicimi funkcemi (Subscript[G, 1]),...,(Subscript[G, n]), pak  ma nahodna velicina (Subscript[X, 1])+...+(Subscript[X, n]) vytvorujici funkci (Subscript[G, 1])(s)×...×(Subscript[G, n])(s). (Dukaz: Stredni hodnota soucinu nezavislych nahodnych velicin je rovna soucinu strednich hodnot.)

Cviceni

Necht (Subscript[X, 1]),...,(Subscript[X, 5]) jsou nezavisle n. v. s geometrickymi rozdelenimi s parametry p=1/6,2/6,...,5/6. Spocitejte vytvorujici funkci pravdepodobnosti n. v. (Subscript[S, 5])=(Subscript[X, 1])+...(Subscript[X, 5]) a pocitejte P((Subscript[S, 5])=150).

⚠ Pred spustenim vypoctu si nezapomente notebook ulozit!

Reseni

tableofgeneratingfunctions=Table[Sum[(s^k)PDF[GeometricDistribution[i/6],k],{k,0,∞}],{i,5}]

Vytvorujici funkce (Subscript[S, 5]) (neni zapotrebi odlozene prirazeni)

generatingfunctionS5[s_]=Apply[Times,tableofgeneratingfunctions]

Prvni pokus:

pst1S5[k_]:=(1/(k!))D[generatingfunctionS5[s],{s,k}]/.s->0

Timing[t90100=Table[pst1S5[k],{k,90,100}];]
t90100//N

Problem je v tom, ze vyssi derivace jmenovatele jsou velmi komplikovane. Pomoci si muzeme funkci Apart (rozklad na parcialni zlomky):

(generatingfunctionS5new[s_]=Apply[Times,tableofgeneratingfunctions]//Apart)//TraditionalForm

Porovnejte slozitost:

{D[generatingfunctionS5[s],{s,2}],
D[generatingfunctionS5new[s],{s,2}]}
%//Simplify

pst1S5new[k_]:=(1/(k!))D[generatingfunctionS5new[s],{s,k}]/.s->0

Timing[t90100new=Table[pst1S5new[k],{k,90,100}];]
t90100new//N

t90100-t90100new

Rozlicna cviceni

1.

1. Pomocí nahrazovacího pravidla proveďte substituci (y^(-1 ))za x ve výrazu: (x^2)Exp[-(x^(-2))].

2.

2. Spočtěte Exp[x] pro záporné řešení rovnice (x^3)-2(x^2)+1==0.

3.

3. Graficky znazornete konvergenci/divergenci harmonicke rady.

4.

4. Prevedte (∑_(n=1)^7)(1/(x+n)) na spolecneho jmenovatele.

5.

5. Najdete reseni soustavy rovnic x=y+1, (y^2)=z,x+z=7.

6.

6. Vyreste rovnici sin,x+cos,x=x.

7.

7. Odhadnete, k cemu konverguje rada (4/π)(∑_(n=1,3,5,...)^∞)((sin,n,x)/n) . Graficky znazornete.

MaxMemoryUsed[]

8.

8. Graficky znazornete zavislost casove narocnosti vypoctu inverzni matice k Hilbertove matici na radu matice n uzitim funkce Timing pro 40≤n≤60. Na vodorovne ose uvedte skutecne hodnoty n. Vypocet nekolikrat zopakujte. (Vysledek demonstruje jev Skienou nazvany WolfRAM, RAM znamena Random Access Memory)

MaxMemoryUsed[]

20130403 Cviceni 6

SetDirectory[NotebookDirectory[]]

AlarmClock

Hyperlink[Button["Budík"],"http://demonstrations.wolfram.com/AlarmClock/"]

FAQ

Rozlicna cviceni 1

1. (3 minuty)

1. Pomocí nahrazovacího (substitucniho, lokalniho) pravidla proveďte substituci (y^(-1 ))za x ve výrazu: (x^2)Exp[-(x^(-2))].

(x^2)(E^(-(x^(-2))))/.x->1/y

2. (5 minut)

2. Spočtěte Exp[x] pro záporné řešení rovnice (x^3)-2(x^2)+1==0.

(* vzor, pattern *)Cases[x/.Solve[(x^3)-2(x^2)+1==0,x],_?Negative]//First//Exp
%//N

nebo

(* vzor, pattern vytvoreny z predikatu *)Cases[x/.Solve[(x^3)-2(x^2)+1==0,x],_?(#<0&)]//First//Exp
%//N

nebo

(* predikat *)Select[x/.Solve[(x^3)-2(x^2)+1==0,x],Negative]//First//Exp
%//N

nebo

(* predikat *)Select[x/.Solve[(x^3)-2(x^2)+1==0,x],#<0&]//First//Exp
%//N

Solve resi algebraicke rovnice (nebo rovnice, ktere lze na algebraicke prevest) presne, NSolve numericky:

Solve[(x^5)-(x^3)-2(x^2)+1==0,x]
%//N

Vybereme kladne koreny:

NSolve[(x^6)-(x^3)-2(x^2)+1==0,x]
Cases[x/.NSolve[(x^5)-(x^3)-2(x^2)+1==0,x],_?Positive]//Exp

3. (10 minut)

3. Graficky znazornete konvergenci/divergenci harmonicke rady.

tabulkaharmonickarada=Table[{m,(∑_(n=1)^m)(1./n)},{m,10,10000,10}];
ListPlot[tabulkaharmonickarada]

Je znamo, ze castecne soucty harmonicke rady se priblizne chovaji jako logaritmus. Graficky znazornime.

fitharmonickarada=Fit[tabulkaharmonickarada,{1,Log[x]},x]

(Vyrovnava funkci α+β, ln(x) metodou nejmensich ctvercu: odhady: (OverHat[α])=0.586191, (OverHat[β])  =0.998953.)

p1=ListPlot[tabulkaharmonickarada]
p2=Plot[fitharmonickarada,{x,10,10000},PlotStyle->{Orange,DotDashed,Thick}]
(* Oba grafy spolecne *)Show[p1,p2]

4. (5 minut)

4. Prevedte (∑_(n=1)^7)(1/(x+n)) na spolecneho jmenovatele.

Zname jiz funkci Apart. Pomoci F1 se podivame nakonec navodu a uvidime “See Also”. Podivame se na Help uvedenych funkci.

Apart

soucet7=(∑_(n=1)^7)(1/(x+n))

Together[soucet7]

Funguje take

Factor[soucet7]

5. (10 minut)

5. Najdete reseni soustavy rovnic x=y+1, (y^2)=z,x+z=7.

Clear[x,y,z]
{x,y,z}/.Solve[{x==y+1,(y^2)==z,x+z==7},{x,y,z}]
{x,y,z}/.Solve[x==y+1&&(y^2)==z&&x+z==7,{x,y,z}]

6. (5 minut)

6. Vyreste rovnici sin, x+cos, x=x.

Plot[Sin[x]+Cos[x]-x,{x,-2π,2π},PlotLabel->Sin[x]+Cos[x]-x]

FindRoot[Sin[x]+Cos[x]==x,{x,-100}]

Poznamenejme, ze Solve a NSolve najde vsechna reseni, FindRoot najde jenom jeden, jeste navic zavisly na nami zvolenem pocatecnim reseni. Ukazeme na nasledujici rovnici:

FindRoot[3 Cos[x]==Log[x],{x,2}]

Reseni pro pocatecni hodnoty 1,...25:

Text@Grid[tabulkaFindRootsCosxMinusLogx=Table[{x0,(x/.FindRoot[3 Cos[x]==Log[x],{x,x0}])},{x0,25}],Frame->{All,False}]

Podivejme se, kolik tech reseni je:

(* Vybereme druhy sloupec tabulky *)vsechnykoreny=tabulkaFindRootsCosxMinusLogx[[All,2]]

18.624716143898166`-18.62471614389822`

(* Vsechny ruzne koreny *)
(* Vidime, ze 18.6247 se vyskytuje dvakrat - priblizna reseni se lisi, jejich rozdil je -5.3290705182007514`*^-14 *)
Union[vsechnykoreny]

Za nerozlisitelne budeme povazovat koreny, ktere se lisi na ctvrtem desetinnem miste. Kolik jich je?

Union[vsechnykoreny,SameTest->(Abs[#1-#2]<(10^(-4))&)]
Length[%]

Plot[3 Cos[x]-Log[x],{x,1,30},AxesOrigin->{0,0}]

Je videt, ze jsme nasli vsechna reseni.

7. (8 minut)

7. Odhadnete, k cemu konverguje rada (4/π)(∑_( n=1,3,5,...)^( ∞))((sin(,n ,x))/n)  na intervalu (-π,π). Graficky znazornete.

7. (Ne) Reseni 1(zatim zustava zahadou)

{Boole[OddQ[2]],Boole[OddQ[3]]}

Clear[soucet];
soucet=(4/π)(∑_(n=1)^∞)((Sin[n x])/n)Boole[OddQ[n]]

Toto funguje

(∑_(n=1)^7)((Sin[n x])/n)Boole[OddQ[n]]

ale toto nefunguje

(∑_(n=1)^m)((Sin[n x])/n)Boole[OddQ[n]]

7.  Reseni 2

Jde vlastne o Fourierovu radu funkce sign(x) na intervalu (-π,π).

Abychom si usetrili praci se substituci indexu uzijeme funkci

FindSequenceFunction[{1,3,5,7},n]

Simplify[(4/π)(∑_(n=1)^∞)((Sin[(2n-1)x])/(2n-1)),Assumptions->(-π<x<π)]
Plot[%,{x,-π,π},PlotStyle->{Red,Thick}]

MaxMemoryUsed[]

20130410 Cviceni 7

FinancialData["USD/EUR"]

SetDirectory[NotebookDirectory[]]

AlarmClock

Hyperlink[Button["Budík"],"http://demonstrations.wolfram.com/AlarmClock/"]

FAQ

Rozlicna cviceni 2

8. Zadani

8. Graficky znazornete zavislost casove narocnosti vypoctu inverzni matice k Hilbertove matici na radu matice n uzitim funkce Timing pro 40≤n≤60. Na vodorovne ose uvedte skutecne hodnoty n. Vypocet nekolikrat zopakujte. (Vysledek demonstruje jev (Skienou^⋆) nazvany WolfRAM, RAM znamena Random Access Memory a nezavisi na nejakem generatoru nahodnych cisel!!!)

⋆Skiena, Steven.: Computational Discrete Mathematics: Combinatorics and Graph Theory in Mathematica. Cambridge University Press, September 2002

8. Reseni

MaxMemoryUsed[]

tableInverseHilbertMatricesTiming[startingdimension_,enddimension_]:=Table[{n,Timing[Inverse[HilbertMatrix[n]];]//First},{n,startingdimension,enddimension}]

tableInverseHilbertMatricesTiming[40,60]
ListPlot[%,Joined->True]

MaxMemoryUsed[]

9. 1. Zadani

Naleznete maxima řádků v dane matici. Uvazujte matici

SeedRandom[13];
(maticecviceni9=RandomVariate[DiscreteUniformDistribution[{0,9}],{3,4}])//MatrixForm

Dále naleznete maxima sloupců.

9. 1. Reseni

Map[Max,maticecviceni9]
Map[Max,Transpose[maticecviceni9]]

9. 2 Zobecneni

Zobecnete na libovolne rozmery v predchozim pripade. Vytvorte obecnou proceduru se dvema argumenty, z nichz prvni bude matice a druhy bude udavat, zda vystupem maji byt maxima radku nebo maxima sloupcu.

9. 2 Reseni

Clear[maximaradkunebosloupcu];maximaradkunebosloupcu[matice_?MatrixQ,radkynebosloupce_]:=Which[ToString[radkynebosloupce]=="radky",Map[Max,matice],
ToString[radkynebosloupce]=="sloupce",Map[Max,matice//Transpose],
True,"Nevyhovuje druhy argument"]

maximaradkunebosloupcu[4,"radky"]
maximaradkunebosloupcu[4,{4}]
maximaradkunebosloupcu[{{4}},{4}]
maximaradkunebosloupcu[{{4}},"sloupce"]
maximaradkunebosloupcu[maticecviceni9,"radky"]
maximaradkunebosloupcu[maticecviceni9,"sloupce"]

10. Zadani

Uvazujte nasledujici data {{(Subscript[x, 1]),(Subscript[y, 1])},...,{(Subscript[x, n]),(Subscript[y, n])}}

SeedRandom[13];
data10={Range[50],Map[RandomVariate[NormalDistribution[2#-(((((#)/7)))^2),0.2#]]&,Range[50]]}//Transpose;

10. 1. Data graficky znazornete.

10. 2. Vyrovnejte data a) primkou, b) parabolou

10. 3. Vysledky znazornete ve spolecnem grafu spolu s originalnimi daty.

10. Reseni

plot1=ListPlot[data10,PlotStyle->PointSize[0.02]]

FindFit[data10,a+b x,{a,b},x]
a+b x/.FindFit[data10,a+b x,{a,b},x]

plot2=Plot[(a+b x/.FindFit[data10,a+b x,{a,b},x])//Evaluate,{x,1,50},PlotStyle->{Thick,Orange}]

(Bez Evaluate bychom museli jako prvni argument Plot primo funkci  12.0548 +0.810952 x. Viz nasledujici reseni.)

Alternativne

ffit=Fit[data10,{1, x},x]
Plot[ffit,{x,1,50},PlotStyle->{Thick,Orange}]

Pro parabolu

FindFit[data10,a+b x+c (x^2),{a,b,c},x]

plot3=Plot[(a+b x+c (x^2)/.FindFit[data10,a+b x+c (x^2),{a,b,c},x])//Evaluate,{x,1,50},PlotStyle->{Thick,Green}]

Show[plot1,plot2,plot3]

Data byla zamerne zatizena chybami. Zobrazovana funkce je

(fskutecna=(2#-(((((#)/7)))^2))&[x])//N
plot4=Plot[fskutecna,{x,1,50},PlotStyle->{Thick,Blue}];
(* vsechny grafy pohromade *)
Show[plot1,plot2,plot3,plot4]

Nebo take

Show[plot1,Plot[{(a+b x/.FindFit[data10,a+b x,{a,b},x]),(a+b x+c (x^2)/.FindFit[data10,a+b x+c (x^2),{a,b,c},x]),fskutecna}//Evaluate,{x,1,50},PlotStyle->Thick,PlotLegends->"Expressions"]]

20130417 Cviceni 8

FinancialData["USD/EUR"]

Directory[]
SetDirectory[NotebookDirectory[]]

AlarmClock

Hyperlink[Button["Budík"],"http://demonstrations.wolfram.com/AlarmClock/"]

FAQ

11. Zadani

Vytvorte funkci, ktera z daneho seznamu vybere vsechny prvky lezici v intervalu [a,b].

11. Reseni

vyberprvkySelect1[seznam_?VectorQ,a_?NumberQ,b_?NumberQ]:=Select[seznam,a≤#≤b&]
vyberprvkySelect2[seznam_?VectorQ,a_?NumberQ,b_?NumberQ]:=Select[seznam,Function[🐺,a≤🐺≤b]]
vyberprvkyCases[seznam_?VectorQ,a_?NumberQ,b_?NumberQ]:=Cases[seznam,_?(a≤#≤b&)]
vyberprvkySelect3[seznam_?VectorQ,{a_,b_}]:=Select[seznam,Function[🐺,a≤🐺≤b]]

vyberprvkySelect1[Range[40],20,30]
vyberprvkySelect2[Range[40],20,30]
vyberprvkyCases[Range[40],20,30]
vyberprvkySelect3[Range[40],{20,30}]

12. Zadani

Vytvorte funkci, ktera ma argumentem seznam a vraci hodnotu True, kdyz jsou vsechny prvky seznamu suda cisla a False, je-li tomu naopak. V obou pripadech zjistete pocet sudych cisel v seznamu. Pouzijte funkci And.

12. Reseni

Map[EvenQ,Range[10]]
Count[%,True]

sudacisla[seznam_]:=Print[Apply[And,Map[EvenQ,seznam]],"  Pocet sudych ",Count[Map[EvenQ,seznam],True]]
sudacisla[Range[10]]
sudacisla[2Range[10]]
sudacisla[2Range[10]+1]

13. Zadani

Uvazujme soucasne ceske mince

ceskemince={"kor","dvo","pet","des","dva","pad"};

se zrejmymi denominacemi {1,2,5,10,20,50}. V kapse mate

SeedRandom[13];
kapsa=RandomChoice[ceskemince,45]

1. Zjistete celkovou hodnotu minci, ktere mate v kapse.

2. Zjistete, kolik kusu od jednotlivych denominaci mate (tj. 10 korun, 3 dvoukoruny, ...)

13. Reseni

1. Celkova hodnota

substituce=Thread[ceskemince->(denominace={1,2,5,10,20,50})]
hodnotyminci=kapsa/.substituce
Total[%]

2. Cetnosti

Text@Grid[cetnosti=Tally[kapsa],Frame->All]

Test, zda jsme si v predchozim pripade pocinali spravne:

cetnosti/.substituce
Map[Apply[Times,#]&,%]
Total[%]

Zamyslete se nad tim, proc nefunguje

Map[Times,cetnosti/.substituce]

14. Zadani

Vytvorte proceduru pro vypocet soucasne hodnoty penezniho toku v ekvidistatnich casovych okamzicich. Uvazujte (i) konstatni hodnotici urokovou miru, (ii) promenne (okamzite, spotove) urokove sazby, (iii) promenne urokove sazby v jednotlivych obdobich.

14. Reseni

(1/(((1+i))^(Range[5]-1)))

Clear[pv];
pv[cf_,i_]:=cf.(1/(((1+i))^(Range[Length[cf]]-1)))
pv[{cf0,cf1,cf2,cf3},i]//TraditionalForm
%/.{i->0.05}

20130424 Cviceni 9

FAQ

Konstrukce  _?VectorQ

Vyraz, ktery vyhovuje vzoru a predikat pro nej dava True, predikat ⟶ vzor

Ad Cviceni 7 (Proc nefunguje Boole[OddQ[n]]?)

7. Odhadnete, k cemu konverguje rada (4/π)(∑_( n=1,3,5,...)^( ∞))((sin(,n ,x))/n)  na intervalu (-π,π). Graficky znazornete.

7. (Ne) Reseni 1(zatim zustava zahadou - jiz neplati)

Co jsme ukazovali drive:

{Boole[OddQ[2]],Boole[OddQ[3]]}

Clear[soucet];
soucet=(4/π)(∑_(n=1)^∞)((Sin[n x])/n)Boole[OddQ[n]]

Toto funguje

(∑_(n=1)^7)((Sin[n x])/n)Boole[OddQ[n]]

ale toto nefunguje

(∑_(n=1)^m)((Sin[n x])/n)Boole[OddQ[n]]

Proč?

Pro obecne n OddQ[n] vraci False:

Clear[n];OddQ[n]

takze

Boole[OddQ[n]]

15. Zadani

Zjistete, zda seznam prirozenych cisel {(Subscript[a, 1]),...,(Subscript[a, n])} je permutace.

?Permutation*

15. Reseni

1. zpusob

Map[PermutationListQ,{{4,3,6,1},{3,2,1}}]

2. zpusob

<<Combinatorica` (* Need["Comb ... `"] *)

?PermutationQ

Attributes[PermutationQ]

Map[PermutationQ,{{4,3,6,1},{3,2,1}}]

3. zpusob

Clear[perm2Q];
perm2Q[seznam_?VectorQ]:=(Sort[seznam]==Range[Length[seznam]])
Map[perm2Q,{{4,3,6,1},{3,2,1}}]

Pokusy

list1={4,5,6,1,2,5};
(*
cetnosti1=Tally[list1]
(Max[Map[Last,cetnosti1]]==1)&&
*)
Intersection[Range[Length[list1]],list1]==Range[Length[list1]]

16. Zadani

Generujte m náhodných n-tic z diskrétního rovnoměrného rozdělení s nosičem {1,2,...,n} a zjistěte, jaká je relativní četnost permutací v této množině. Odpovídá to teoretické pravděpodobnosti?

16. Návod

Generování deseti pětic

RandomChoice[Range[5],{10,5}]

Toto se sice k řešení úlohy nepotřebuje, ale je to také užitečná konstrukce:

(* returns a random permutation *)RandomSample[Range[10]]

16. Reseni

tableofmn[m_,n_]:=RandomChoice[Range[n],{m,n}]

tableofmn[10,5]

SeedRandom[1313];
tabulka10000petic=tableofmn[1000000,5];
Count[tabulka10000petic,_?PermutationQ]

((5!)/(5^5))//N

17. Zadani

Simulujte 50 trajektorií ceny akcie řídící se náhodnou procházkou s počáteční cenou 1000 pro 500 obchodních dnů, přičemž denní změna má standardní normální rozdělení N(0,1). Tyto trajektorie znázorněte v jednom grafu.

17. Návod

Návod: nejdříve vytvořte funkci pro jednu náhodnou procházku užitím např. funkce NestList, Accumulate,... . Jedna náhodná procházka = jedna trajektorie. Z těchto jednotlivých trajektorií vytvořte tabulku (tu uchovejte pro další potřebu). Všechny trajektorie znázorněte v jednom grafu.

17. Reseni

jednaprochazka:=NestList[#+RandomVariate[NormalDistribution[]]&,1000,500]
(* Fakticky dostaneme posloupnost i s pocatecni hodnotou, takze je delky 501 *)

SetOptions[ListPlot,Joined->True,PlotStyle->Thick];

Nemáme SeedRandom, po každém spuštění je obrázek jiný. Příklad na třech procházkách:

ListPlot[{jednaprochazka,jednaprochazka,jednaprochazka}]

Volby ListPlot vratit na puvodni:

SetOptions[ListPlot,Joined->False,PlotStyle->Automatic];

Teď už ale pro účely reprodukovatelnosti zasejeme semínko. Nicméně ty náhodné procházky se generují postupně, takže vlastně z jednoho náhodného semínka se generuje 50×500=25000 náhodných čísel! A proto se jednotlivé náhodné procházky liší!

SeedRandom[1313];ListPlot[tabulka50nahodnychprochazekdelky500=Table[jednaprochazka,{50}]]

18. Zadani

Z 50 trajektorií z předchozího příkladu (máte je uchované v tabulce) vytvořte graficky
(i) obálku, tj. křivku maxim a křivku minim ze všech trajektorií,
(ii) průměrnou trajektorii, tj. v každém časovém okamžiku vypočítáte průměr ze všech 50 pozorování (použijete funkci Mean),
(iii) mediánovou trajektorii; místo funkce Mean se použije funkce Median.
Vše znázorněte v jednom grafu.
(iv) graficky znázorněte také volatilitu jako funkci času (pro každý časový okamžik se vypočítá směrodatná odchylka z odpovídajích 50 hodnot.

18. Reseni

(i)

Dimensions[tabulka50nahodnychprochazekdelky500]

My potřebujeme řezy v časových okamžicích:

tabulka501obchodnichdnu=Transpose[tabulka50nahodnychprochazekdelky500];
Dimensions[%]

Obálka maxima:

ListPlot[Map[Max,tabulka501obchodnichdnu]]

ListPlot[{Map[Max,tabulka501obchodnichdnu],Map[Min,tabulka501obchodnichdnu]}]

Kompaktní konstrukce:

ListPlot[Map[Map[#,tabulka501obchodnichdnu]&,{Max,Min,Mean,Median}]]

Volatilita ( = směrodatná odchylka)

ListPlot[Map[StandardDeviation,tabulka501obchodnichdnu]]