My name is Abdul Jawad, and I am a Ph.D. candidate at the
University of California Santa Cruz . My research involves modeling human
driving behavior to assist simulation-based autonomous vehicle testing. My goal is to accelerate the development and deployment of autonomous vehicles by improving their
testing process through advanced simulation techniques. In addition to
publishing at academic conferences, I have contributed to
open-source projects
related to AV simulation tool development. The tools are primarily developed in Python and C++. I have been using Unreal and Unity game engines for seven years due to my exposure to game development.
I am primarily interested in modeling behaviors for simulation-based AV testing for my dissertation project.
I received my BSc in Computer Science and Engineering from Bangladesh. I co-founded a game studio and published four Andriod games
during my undergraduate. In 2018, I moved to the University of California, Santa Cruz, to pursue a PhD in the Computational
Media Department.
Here's what you'll find on this website. Under the
Research tab,
I maintain a list of my publications and research projects.
On the
Projects tab, you can find a bunch of my other projects.
Please find my contact information below this page, as well as other links.
To see all my public project repositories, check out my
GitHub profile.


SmashBoard
Promises, promises, not only words...
Talk to me
Status:
Now playing on satellite TV
I'm not sure boss
Channel #
Location
Latitude:
Longitude:


Or Sign in with social platforms
Or Sign up with social platforms
New to the website? Join now and level up your skills.
Already a member? Log in now and continue learning.


Awwwww...Don't Cry :(
It's just a 404 Error!
What you’re looking for may have been misplaced in Long Term Memory.
Go Home


|
phylopomp
Phylodynamics for POMPs
|
#include "internal.h"
Go to the source code of this file.
Macros | |
| #define | MAKEFN(X, TYPE) |
| #define | REVIVEFN(X, TYPE) |
| #define | RUNFN(X, TYPE) |
| #define | GENEALFN(X, TYPE) |
| #define | YAMLFN(X, TYPE) |
| #define | GENERICS(X, TYPE) |
Functions | |
| template<class TYPE > | |
| SEXP | ndeme (TYPE &X) |
| template<class TYPE > | |
| SEXP | nsample (TYPE &X) |
| template<class TYPE > | |
| SEXP | timezero (TYPE &X) |
| template<class TYPE > | |
| SEXP | time (TYPE &X) |
| template<class TYPE > | |
| SEXP | serial (const TYPE &X) |
| binary serialization More... | |
| template<class TYPE > | |
| SEXP | yaml (const TYPE &X) |
| human/machine readable output More... | |
| template<class TYPE > | |
| SEXP | describe (const TYPE &X) |
| human readable output More... | |
| template<class TYPE > | |
| SEXP | structure (const TYPE &X) |
| structure in R list format More... | |
| template<class TYPE > | |
| SEXP | newick (const TYPE &X) |
| tree in newick format More... | |
| template<class TYPE > | |
| SEXP | lineage_count (const TYPE &G) |
| number of lineages through time More... | |
| template<class TYPE > | |
| SEXP | make (SEXP Params, SEXP IVPs, SEXP T0) |
| initialization More... | |
| template<class TYPE > | |
| SEXP | revive (SEXP State, SEXP Params) |
| refresh parameters More... | |
| template<class TYPE > | |
| SEXP | run (SEXP State, SEXP Tout) |
| run simulations More... | |
| template<class TYPE > | |
| SEXP | genealogy (SEXP State) |
| extract the bare genealogy More... | |
| #define GENEALFN | ( | X, | |
| TYPE | |||
| ) |
Definition at line 135 of file generics.h.
| #define GENERICS | ( | X, | |
| TYPE | |||
| ) |
Definition at line 143 of file generics.h.
| #define MAKEFN | ( | X, | |
| TYPE | |||
| ) |
Definition at line 123 of file generics.h.
| #define REVIVEFN | ( | X, | |
| TYPE | |||
| ) |
Definition at line 127 of file generics.h.
| #define RUNFN | ( | X, | |
| TYPE | |||
| ) |
Definition at line 131 of file generics.h.
| #define YAMLFN | ( | X, | |
| TYPE | |||
| ) |
Definition at line 139 of file generics.h.
| SEXP describe | ( | const TYPE & | X | ) |
human readable output
Definition at line 44 of file generics.h.
| SEXP genealogy | ( | SEXP | State | ) |
extract the bare genealogy
Definition at line 114 of file generics.h.
| SEXP lineage_count | ( | const TYPE & | G | ) |
number of lineages through time
Definition at line 62 of file generics.h.
| SEXP make | ( | SEXP | Params, |
| SEXP | IVPs, | ||
| SEXP | T0 | ||
| ) |
| SEXP ndeme | ( | TYPE & | X | ) |
Definition at line 7 of file generics.h.
| SEXP newick | ( | const TYPE & | X | ) |
tree in newick format
Definition at line 56 of file generics.h.
| SEXP nsample | ( | TYPE & | X | ) |
Definition at line 12 of file generics.h.
| SEXP revive | ( | SEXP | State, |
| SEXP | Params | ||
| ) |
refresh parameters
Definition at line 87 of file generics.h.
| SEXP run | ( | SEXP | State, |
| SEXP | Tout | ||
| ) |
| SEXP serial | ( | const TYPE & | X | ) |
binary serialization
Definition at line 28 of file generics.h.
| SEXP structure | ( | const TYPE & | X | ) |
structure in R list format
Definition at line 50 of file generics.h.
| SEXP time | ( | TYPE & | X | ) |
| SEXP timezero | ( | TYPE & | X | ) |
Definition at line 17 of file generics.h.
| SEXP yaml | ( | const TYPE & | X | ) |
human/machine readable output
Definition at line 38 of file generics.h.


[enable JavaScript to render MarkDown]:
# Features of Scrollback ####**Infinite scroll** Scroll up to see what people were talking about in the last hour - or last week. ####**Minimize ticker** Gently introduce casual visitors to your community with conversation that is never out of sight. ####**Subtle notifications** Notifications and timestamps appear when needed, then fade away to keep the focus on the conversation. ####**User colors** Quickly identify who's talking to whom, with colors that represent people consistently. ####**Timeline** See the peaks and troughs of conversation in a colorful visualization you can scroll through. ####**User history** Hover over a nickname to see that person's activity highlighted on the timeline. ####**IRC compatible** If an active IRC community exists, we gateway to it instead of competing with it. ####**Seamless integration** Let Scrollback match the fonts and styles of your website automatically, or customize it extensively with CSS.


Congratulations, this page is being served up by a simple HTTP server running from your Nova compute instance!
Common OpenStack Deployments is brought to you by Elizabeth K. Joseph and Matt Fischer.
Welcome to the...
Congratulations, this page is being served up by a simple HTTP server running from your Nova compute instance!
Common OpenStack Deployments is brought to you by Elizabeth K. Joseph and Matt Fischer.


Enter a string to get its FNV hash:
Not associated with or endorsed by Electronic Arts or Maxis.
Created by Kade. Contact via Discord or email. Source code on GitHub.


Track Selector
So you want to learn how to program? Take our short survey by cllicking the button below and we'll help you choose the correct track.
C#
C# is most popular among bigger established businesses, often for building internal software. Because it's been around for a long time and has the backing of Microsoft, it is one of the most in-demand languages in the job market. C# has also been going through a bit of a rebirth lately, with Microsoft open sourcing the language and surrounding platform, porting it to run on Mac and Linux, and incorporating many of the best features of other languages. If you like the idea of working for a larger company on business software, C# is a great choice.
Ruby
Ruby is a favorite language of developers building interactive web applications. If an app involves users creating accounts, entering information, and interacting with dynamic content, there's a good chance it is built with Ruby. Ruby became popular because the Rails framework, which is written with Ruby, simplified many of the common tasks associated with building web applications. It's most popular with startups and smaller companies who are looking to build their product quickly.
PHP
PHP is by far the most popular backend language today, with 80 percent of websites utilizing it 'server-side'. It is perhaps best known for it's use in content management systems like Wordpress, Drupal, and Joomla. But the versatility of the language and the frameworks it powers make employment options numerous and diverse. If you're keen to work for a fast paced agency that builds websites for lots of clients, or maintain the security and stability of a huge complex of government websites, or if you like the idea of building out small sites for brands, businesses, and organizations - In any of these cases, PHP would be a great way to go.


This is a guide for making IDPF Epub ebooks/publications. It is mostly an annotated example: this document itself in Epub form.
Not all details/variations are mentioned, but enough to obviate need of the specifications for normal use. And it is for making entirely conformant publications.
Included also is a description of optional extra styling for a particular reader (but still completely conformant).
You need an understanding of and ability to make XHTML/CSS and XML documents.
‘Epub’ is a standard from the International Digital Publishing Forum. It is an arrangement of several other standards (mainly: XHTML, CSS, XML, NCX, DCMI). There are three parts, addressing: content, package metadata, and archive (OPS, OPF, and OCF). It is powerful, straightforward, and non-proprietary.
‘ADE’ is one of the first readers for Epub publications. It is very conformant with the standard. It can use an optional proprietary publication component: an extra stylesheet to adjust text-column appearance. (That is allowed by the standard.)
This guide was written using ADE version 1.0.467 .


--- layout: index title: Home permalink: / ---
This site contains an integrated collection of tools to support future water resources planning and management. This includes advances in hydrometeorological datasets, climate downscaling, hydrologic modeling, climate impact assessements, and streamflow prediction, packaged as open-source software, freely available datasets, and guidance on appropriate use. Collectively these tools provide a foundation for understanding and adapting to future change. This website is dedicated to explaining each tool and how they fit together, documenting for both the science and applications communities resources that are currently available and major activities currently underway.
The mission of the National Center for Atmospheric Research's (NCAR) Research Applications Laboratory is to conduct directed research that contributes to the fundamental understanding of the atmosphere and related physical, biological, and social systems; to support, enhance, and extend the capabilities of the scientific community, and to develop and transfer knowledge and technology for the betterment of life on Earth


BATMAN A film review by Randy Parker Copyright 1996 Randy Parker
RATING: ***1/2 (out of ****)
(Review written in 1989)
What surprises me most about the long-awaited BATMAN is that the movie centers not on the Caped Crusader but on the Joker. The Dark Knight actually has less screen time than his homicidal arch nemesis! I can hear Batfans everywhere protesting: "First they had the nerve to cast Michael Keaton as Batman, and now they think they can get away with making him second banana to the Joker!!" And, indeed, this situation would have been altogether inappropriate and intolerable were it not for the fact that someone at Warner Brothers had the savvy to cast Jack Nicholson in the role he was born to play.
Nicholson deserves the top billing he is receiving in the posters and opening credits since his astounding, Oscar-worthy performance as the Joker is the highlight of BATMAN. Nicholson shows no inhibition or shame as he pulls out all the stops to create a Joker that makes Cesar Romero look like Bozo the Clown (Romero played the Joker in the '60s TV series). On one hand, Nicholson's Joker is extremely funny. His murderous escapades and manic outbursts (which include defacing paintings in an art museum and invading the television airwaves with hilarious pirate broadcasts) bring humor and laughter into an otherwise dark film. But on the other hand, his outwardly jovial disposition is only a mask for his intense contempt for all of mankind. The Joker's cruelty and demented, diabolical nature are truly frightening and sinister.
But enough about Nicholson as the Joker: I'm sure most people are wondering about Michael Keaton as Batman. Like many, I was outraged when Warner Brothers announced that Keaton had been cast as the Caped Crusader. But Batfans need not worry becauseholy unconventional casting!Keaton delivers. He is convincing as the tormented millionaire, Bruce Wayne, and beneath his hard, black body armor, he is mysterious and formidable as the Batman. Even though the screenplay only halfheartedly attempts to develop the character, Keaton manages to give hints of the inner turmoil and emotional rage that drive Wayne to don his Bat gear and battle the low life scum of Gotham City.
Visually, BATMAN is a knockoutnothing short of a triumph for director Tim Burton and his production team. The movie successfully captures the film noir look and bleak, gloomy atmosphere of such Batman comic books as "The Dark Knight Returns" and "The Killing Joke." Burton's film is a feast for the eyes; the $35 million budget is definitely on the screen. The sets, costumes, and make-up are stunning and worthy of Oscar nominations.  His grimacing smile, green hair, and purple suits over orange shirts are colorful and outrageous. Keaton's wardrobe is equally appropriate. Batman's cape, mask, and body armor make him look threatening, frightening, and invinciblebringing to mind not Adam West's campy costume in the TV series but Peter Weir's outfit in ROBOCOP.
Production designer Anton Furst's sets are a wonder. Gotham City is dark and foreboding, almost as if it were an expressionist vision of New York. The Batcave is cold, bare, and desolate. The world that Burton and Furst have created in BATMAN is a melancholy nightmare of violence and corruption. The movie's gritty BLADERUNNER/ROBOCOP look is a joy to behold.
The plot, on the other hand, is not as rousing. The screenplay is not up to snuff in comparison to the stylized sets and superb acting. It takes forever for what little story there is to finally kick into gear. BATMAN tries to walk a thin line between action and character psychology, but it fails to give us enough of either one. As I mentioned earlier, the screenplay gives us plenty of the Joker, but it never fully establishes Batman's psychological motives and underpinnings. In terms of action, BATMAN is not as thrilling as you might expect. Burton seems more concerned with the visual and dramatic elements of the story.
While BATMAN may not completely absorb and excite you, it will likely impress you with its success at bringing to life the style and spirit of the Batman comic books. Thanks to its distinctive look, the movielike a good comic booktransports us to another world.
~~~~~~~~~~ Randy Parker http://www.shoestring.org --> "1,000 Voices In The Dark" --> Weekly Movie Reviews & Interviews
The review above was posted to the
rec.arts.movies.reviews newsgroup (de.rec.film.kritiken for German reviews).
The Internet Movie Database accepts no responsibility for the contents of the
review and has no editorial control. Unless stated otherwise, the copyright
belongs to the author.
Please direct comments/criticisms of the review to relevant newsgroups.
Broken URLs inthe reviews are the responsibility of the author.
The formatting of the review is likely to differ from the original due
to ASCII to HTML conversion.
Related links: index of all rec.arts.movies.reviews reviews


${summary_title}
Board: ${board_type}
Summary
${table_summary_results}
Full Results
#
Result
Name
Description
Time
${table_results}
Total test time: ${total_test_time}.
See
full test output
or
see parameters and retest
.
Station: ${location} ${station}.
Questions? Email
boardfarm-admins
.


Request Payment of Digital Assets
To address:
To tag (required for some chains):
Asset:
Amount:
Description (optional):
UniqueTxId (optional):
Request
Get account to send assets to user
getAccount


Blockchain for Business
WITH IBM BLOCKCHAIN
TRACK & TRACE
Home
Manufacturer
Distributor
Dealer
Hospital
Dashboard
Could not update data
SUCESSFULLY UPDATED PRODUCT STATUS
No product with specified uuid
Dashboard
UUID
Received
Received
Received
Received


Welcome to FirebaseRTC!
perm_camera_mic
Open camera & microphone
group_add
Create room
group
Join room
close
Hangup
Join room
Enter ID for room to join:
Room ID


Module dynamixel-bus
Library for Dynamixel protocol.
This library allows to manipulate devices that use Dynamixel
protocol, such as AX-12 robotic servo motors.
When available, a dynamixel bus will be represented by a Device
object in toribio.devices table. The device will be named (as an
example), "dynamixel:/dev/ttyUSB0".
| motors |
Motors connected to the bus. |
| events |
Signals emitted by this device. |
| conf |
Configuration Table. |
| name |
Name of the device. |
| module |
Module name (in this case, 'dynamixel'). |
| filename |
Device file of the bus. |
| discovery_finished |
The bus scannning status. |
| task |
Task that will emit signals associated to this device. |
Functions
-
init (conf)
-
Initialize and start the module.
This is called automatically by toribio if the load attribute for the module in the configuration file is set to
true.
Parameters:
- conf
the configuration table (see conf).
-
get_broadcaster ()
-
Get a broadcasting Motor object.
All commands sent to this motor will be broadcasted
to all motors.
Returns:
A Motor object.
-
get_motor (id)
-
Get a Motor object.
Parameters:
- id
The numeric ID of the motor
Returns:
A Motor object, or nil if not such ID found.
-
get_sync_motor (...)
-
Get a Sync-motor object.
A sync-motor allows to control several actuators with a single command.
The commands will be applied to all actuators it represents.
The "get" methods are not available.
Parameters:
- ...
A set of motor Device objects or numeric IDs
Returns:
a sync_motor object
-
has_errors (code)
-
Decodes dynamixel error codes.
Parameters:
- code
A dynamixel error code as returned by the different motor methods.
Returns:
A set of error strings. The possible error strings are 'ERROR_INPUT_VOLTAGE', 'ERROR_ANGLE_LIMIT',
'ERROR_OVERHEATING', 'ERROR_RANGE', 'ERROR_CHECKSUM', 'ERROR_OVERLOAD', 'ERROR_INSTRUCTION'.
The table is double-indexed by entry number to allow array-like traversal.
Usage:
local errorset = dynamixel.has_errors(0x10+0x08)
print ("has errors:" #errorset>0 )
print ("has range error:", errorset['ERROR_RANGE']==true )
Tables
-
motors
-
Motors connected to the bus.
The keys are device numbers, the values are Motor objects.
-
events
-
Signals emitted by this device.
Fields:
- discovery_end
The discovery of attached motors has finished.
-
conf
-
Configuration Table.
When the start is done automatically (trough configuration),
this table points to the modules section in the global configuration table loaded from the configuration file.
Fields:
- load
whether toribio should start this module automatically at startup.
- filename
the device file for the serial bus (defaults to '/dev/ttyUSB0').
- serialtimeout
the timeout when waiting for a response from an actuator (defaults to 0.05s)
- stty_flags
parameters to pass to stty for configuring the serial bus.
- serialspeed
serial bus speed in bps (defaults to 1000000).
Fields
-
name
-
Name of the device.
Of the form 'dynamixel:/dev/ttyUSB0'
-
module
-
Module name (in this case, 'dynamixel').
-
filename
-
Device file of the bus.
For example, '/dev/ttyUSB0'
-
discovery_finished
-
The bus scannning status.
This flag will be set to true once the bus has been scanned and all motors discovered.
-
task
-
Task that will emit signals associated to this device.


-
AmazingOutlet
Clean & Elegant design with a modern style. This template includes
all you need for a fashion & accessories store
Shop Now !
-
Women Clothing
Clean & Elegant design with a modern style. This template includes
all you need for a fashion & accessories store
Shop Now !


fswebcam filename- Take a picture with custom resolution:
fswebcam -r widthxheight filename- Take a picture from selected device(Default is `/dev/video0`):
fswebcam -d device filename- Take a picture with timestamp(timestamp string is formatted by strftime):
fswebcam --timestamp timestamp filename


@import model.Institution ; var prevmd = ""; @(item: model.Item, institution: Institution, community: model.Community) @main(institution, community) {


+obj
+obj
is a base class used to turn a DB object into a first
class GUI item. See +Obj, +ObjView and +ObjVal.
+Obj
+Obj
supplies suggestions from the database, using
dbHint
.
(gui '(+E/R +Obj +TextField) '(pa : home obj) '(nm +Man) 30)
Here, the +TextField
is connected to the 'pa'- property of the
current form object. Since we are using +Obj
in combination
with +TextField
, the database will suggest names from the 'nm'-
property of the class +Man
. See also +NumField, +DateField, dbHint, +hint, +obj.
+ObjView
+Obj
, but not editable (view-only). Therefore it
is free to display an arbitrary attribute of an object, or even the result
of an arbitrary formatting expression, while still allowing to click on it
and jump to the object's GUI form, as determined by the url>- method (or
a custom method supplied as an optional argument to +ObjView
.
Such an optional method is accepted by all subclasses of +obj
).
A common use of +ObjView
is when presenting data in a
chart
. Given the following 'put'- function:
'((This) (list This This (: pa) (: ma) (: mate) (: job) (: dat)))
The data can be displayed with the following table:
(<table> 'chart NIL
'((btn) (NIL "Name") (NIL "Father") (NIL "Mother") (NIL "Partner") (NIL "Occupation") (NIL "born"))
(do (cho)
(<row> (alternating)
(gui 1 '(+DstButton) Dst)
(gui 2 '(+ObjView +TextField) '(: nm))
(gui 3 '(+ObjView +TextField) '(: nm))
(gui 4 '(+ObjView +TextField) '(: nm))
(gui 5 '(+ObjView +TextField) '(: nm))
(gui 6 '(+TextField))
(gui 7 '(+DateField)) ) ) )
Columns 2, 3 and 4 use +ObjView
to display the name of the 'pa',
'ma' and 'mate' relations of the current object respectively. See also +Obj, +obj.
+ObjVal
+Obj
, but automatically creates a new object if one
does not exist and an attribute value is entered into a field. See also +Obj, +obj.
+OnClick
(gui '(+OnClick +Button) "alert('OK')" "Alert")
The PicoLisp wiki includes a good article that
explains +OnClick in depth. See also +Button
and gui.
(okButton 'Exe)


Stochastic optimization I: from Monte-Carlo methods to Gibbs sampling
| Date: | Feb 15, 2019 |
|---|---|
| Last Updated: | Feb 20, 2019 |
| Categories: | |
| Tags: |
Contents
Introduction
This is a series of inspection on stochastic methods. Traditionally, an optimization problem could be solved by gradient descent methods, greedy algorithm and stochastic methods. For example, Levenberg–Marquardt algorithm is a gradient descent based methods. Another example is OMP which is used to find a local minimum of L0 penalized problem. Since the L0 norm is totally indifferentiable, we could not calculate its gradient. Therefore, such kind of problem is generally more difficult than those differentiable problems. However, in this series, we would only talk about the stochastic methods. When applying this method to solve optimization, we usually need to construct a probability which shows how possible the solution lies in a specific location.
In this series, we would starter from a beginner’s view on stochastic optimization, and end on some state-of-art researches. A person who read this article may need to have fundamental experiences on such issues:
Gradient based optimization.
Probability theory on elementary level.
Stochastic process on elementary level.
Note that in many cases, a method should not be classified into a particular category, because it may integrate some different ideas from different ares. For example, to solve a Bayesian analysis problem, we could apply the gradient descent method on the likelihood function instead of using random search. Generally, if you know more about optimization, you may have a stronger basis for learning stochastic methods. Actually, many newest topics including deep learning, entropy analysis and uncertainly quantification are related to stochastic methods. In this series, we may talk about those topics.
Certainly, in this topic, we would have some learning materials including lecture notes, slides, blogs and papers. To understand some concepts better, we may provide some codes. In the following notes, we are mainly focus on some concepts that are not explained well in our reading materials. Some topics mentioned in those materials but not very important would be ignored.
This topic
In this topic, we would discuss Monte-Carlo method, Metropolis-Hastings algorithm and Gibbs sampling. The learning materials could be found here:
Monte-Carlo method
What is Monte-Carlo method
Check here to see Monte Carlo method in Wikipedia.
Monte-Carlo method is a general experiment for estimating a integration. For example, for any function $f(x)$, we could calculate the integration by
where $x_i \sim U(a,b)$, which means if we generate a lot of samples that follow the uniform distribution, then we could estimate the integration statistically. Since $f(\cdot)$ is an arbitrary function, we know that a uniform sampling could be used to estimate arbitrary integration. This idea is easy to be extended to the case of vector calculation.
Furthermore, to get a more general form, we could rewrite $\eqref{fml:mc:intro}$ as
where $x_i \sim p(\cdot)$. $p(\cdot)$ is an arbitrary distribution, which means we could use an arbitrary distribution to generate a lot of samples, and use those samples to estimate an arbitrary integration.
| An example of estimating $\pi / 4$ | An example of Buffon’s needle test |
|---|---|
Check here to see Buffon's needle problem in Wikipedia.
Here are two examples for estimating $pi$ by Monte-Carlo method. In the first example, we count the points that satisfy $\lVert (x,y) \rVert_2 < 1$ and use the ratio of these points to estimate $\frac{\pi}{4}$. In the second example, we throw needles randomly and count the ratio of needles that lie across the parallel lines. The ratio could be used to estimate $\pi$. In both of these examples, we do not calculate any integration, and we do not need to measure any area. However, we could use a lot of experiments to approximate the unknown integration, $\pi$.
Unbiased estimator and error analysis
Check here to see Bias of an estimator in Wikipedia.
Unbiased estimator for the expectation
Remark that $\eqref{fml:mc:general}$ is called “Expectation over distribution $p(\cdot)$“. Note that if we denote that
where $\mu$ is the expectation and $\overline{f(x)}$ is the estimator for $\mu$. Because for any $i$, we have $\mathbb{E}_p[f(x_i)] = \mu$, we have
So $\overline{f(x)}$ is the unbiased estimator for the expectation $\mu$. Note that although $\overline{f(x)}$ is a random variable with an average of $\mu$, it is not equivalent to the real expectation. However, this derivation shows that why we could use the average of stochastic samples to replace the original integration in $\eqref{fml:mc:general}$.
Unbiased estimator for the variance
Consider the standard deviation $\sigma$, we know that $\sigma:=\sqrt{\mathbb{D}(\cdot)}$, where $\mathbb{D}$ is the variance. Then if we have a distribution $p(\cdot)$, we would have
First, we have already known that for any function $f$, $\mathbb{D}_p[f(x)] = \mathbb{E}_p \left[ f^2(x) \right] - \mathbb{E}^2_p \left[ f(x) \right]$. Because for any $i \neq j$, $\mathbb{E}_p \left[ f^2(x_i) \right] = \mathbb{E}_p \left[ f^2(x) \right]$ and $\mathbb{E}_p \left[ f(x_i)f(x_j) \right] = \mathbb{E}^2_p \left[ f(x) \right]$, we could derive that
Hence $s^2$ is an unbiased estimator for the real variance $\sigma^2$, because the expectation of $s^2$ is exactly $\sigma^2$.
Check here to see Unbiased sample standard deviation in Wikipedia.
A following question is: “is $s$ also the unbiased estimator for $\sigma$?” Centainly not. Assume that we have already known $\mathbb{E}(x^2) = a$, generally $\mathbb{E}(x) \neq \sqrt{a}$, unless $x$ is a constant. In fact $s$ is not a constant, so $s$ is not the unbiased estimator for $\sigma$.
The second question is: “Is it possible for us to calculate the definite value of the unbiased estimator for the standard deviation?” In some cases we could do that. Since this technique requires us to use Taylor expansion for $\sqrt{\cdot}$, it is very difficult for most cases. The reference from Wikipedia here gives some examples of the unbiased estimator for standard deviation, but here we would not discuss about this topic further.
Unbiased estimator for the error of expectation
Check here to see Standard error in Wikipedia.
Now we would introduce a new concept: standard error. It is important to know that STANDARD ERROR is quite different from STANDARD DEVIATION. A deviation is determined by the specific distribution $p(\cdot)$ and the function $f(\cdot)$. Once we determine $p(\cdot),~f(\cdot)$, the $\sigma$ would be a constant value. Although we could use the unbiased estimator $s^2$ to estimate $\sigma^2$, even we generate infinite samples, the estimated variance would only converge to $\sigma^2$ (not 0). However, the standard error means how much the variance of the estimated expectation is.
Recall that we use $\overline{f(x)} := \frac{1}{N} \sum_{i=1}^N f(x_i)$ to make an unbiased estimation, and we know that such an estimator would converge to the true expectation exactly, i.e. $\mathbb{E}_p\left[\overline{f(x)}\right] = \mu$ (proved in $\eqref{fml:mc:unb-exp}$). But we know that $\overline{f(x)}$ is also a random variable, so it would has the variance (deviation). Hence we denote the standard deviation of the estimator for the expectation as the standard error, i.e.
According to the definition of the variance in $\eqref{fml:mc:def-var}$, for any $i$, we have $\mathbb{E}_p \left[ \left( f(x_i) - \mu \right)^2 \right] = \sigma^2$, where $x_i \sim p(\cdot)$. Then we have
In the above derivation, $\varepsilon^2 = \frac{1}{N} \sigma^2$ discloses that standard error which is a constant is exactly $\frac{\sigma}{\sqrt{N}}$, where the standard deviation is also a constant. Hence, if $\mathbb{E}_p \left[ s^2 \right] = \sigma^2$ shows that $s^2$ is the unbiased estimator for the variance $\sigma^2$, then we know the unbiased estimator for $\varepsilon^2$ is
It is easy to know that
To be emphasized, $\mathrm{SE}^2$ is an unbiased estimator for $\varepsilon^2$, however, $\mathrm{SE}$ is not the unbiased estimator for $\varepsilon$. Because we are borrowing the derivation of the estimator for the variance. Certainly, if we could calculate the unbiased estimator for $\sigma$ (although it is difficult to find such an estimator), then it is easy to learn that the unbiased estimator for $\varepsilon$ is $\frac{1}{\sqrt{N}}$ of the estimator for $\sigma$.
Importance sampling
Defining a set (area) that $\mathcal{Q} = \{ x | q(x) > 0 \}$, we know it has a subset $\mathcal{X} = \{ x | f(x)q(x) \neq 0 \}$. Consider another set $\mathcal{P} = \{ x | p(x) > 0 \}$, then we know that $\mathcal{P}$ and $\mathcal{Q}$ are the domains of two different distributions $p(\cdot),~q(\cdot)$ respectively. As long as we ensure that $\mathcal{X} \in \mathcal{P}$, we could calculate the expectation $\mathbb{E}_q [f(x)]$ by
| General case | Simplified case |
|---|---|
The above two figures show the relationship among three domains. In the left figure, we show a general case, where $\mathcal{P} \neq \mathcal{Q}$. In some cases, when $q(x) > 0$, $p(x) = 0$, but we have $f(x)=0$ in such cases. So we could use a distribution like $p(\cdot)$ with a domain $\mathcal{P}$ to calculate the expectation.
Such a relationship is not convenient in many cases. To simplify it, we could construct a distribution $p(\cdot)$ which domain satisfies $\mathcal{Q} \subseteq \mathcal{P}$. This case is shown in the right figure.
The technique where we use $\mathbb{E}_p \left[f(x)\frac{q(x)}{p(x)}\right]$ to calculate $\mathbb{E}_q [f(x)]$ is called “importance sampling” which means we use a distribution to calculate Monte-Carlo integration over another distribution.
Unbiased estimator for importance sampling
Since $\mathbb{E}_p \left[f(x)\frac{q(x)}{p(x)}\right] = \mathbb{E}_q [f(x)]$, if we denote $w_i = \frac{q(x_i)}{p(x_i)}$, we could define such unbiased estimators
where $x_i \sim p(\cdot)$.
If we set $f^{\dagger}(x_i) := w_i f(x_i)$, we would know that the three estimators could be proved to be unbiased by the totally same methods discussed before. These three estimators are for expectation, variance and squared standard error respectively.
To be specific, we have known that $\mathbb{E}_p \left[f(x)\frac{q(x)}{p(x)}\right] = \mathbb{E}_q [f(x)]$, which means $\mathbb{E}_p \left[ \overline{wf(x)} \right] = \mathbb{E}_q [f(x)] = \mu_q$. So the unbiased estimator for the expectation is exactly the real expectation over the $q(\cdot)$.
However, If we consider the unbiased estimator for the variance, we would find that
Consider the real variance $\sigma^2_q = \mathbb{E}_q \left[ f^2(x) \right] - \mu^2_q$, we would find that the bias of the proposed variance is
Hence we know that the unbiased estimator for the “variance” here is not the unbiased estimator for the true variance of $f(\cdot)$ over $q(\cdot)$. Certainly, it is really an unbiased estimator, but it converges to another variance.
Since $\exists~C,~\forall~x, \frac{q(x)}{p(x)} < C$, we know that $\widehat{\sigma_q^2}$ would not be extremely larger than $\sigma_q^2$ if we choose a good density function $p(\cdot)$ that is close to $q(\cdot)$. In particular, when $p(\cdot) = q(\cdot)$, $\widehat{\sigma_q^2} = \sigma_q^2$.
There is the similar circumstance for the estimator for standard error, because the $\widehat{\varepsilon^2_q} = \frac{1}{N} \widehat{\sigma_q^2}$.
The above derivations shows that when using importance sampling to calculate the Monte-Carlo integration, we could still preserve the expectation for another probability. However, compared to generating samples over the original probability directly, this method may introduce a different variance (almost always larger). The above unbiased estimators could estimate the variance and standard error correctly, but both of the variance and standard error are different from those of the direct sampling method.
Alternative: self-normalized importance sampling
Check here to see Importance sampling in Wikipedia.
We have already known that $\overline{wf(x)}$ is an unbiased estimator for $\mathbb{E}_q [f(x)]$. Let $f(x) = 1$, then we could define an unbiased estimator
where $\mathbb{E}_p (\overline{w}) = 1$.
Since
We could derive an estimator that
According to $\mathbb{E}\left[ \frac{a}{b} \right] = \frac{\mathbb{E}[a]}{\mathbb{E}[b]} - \frac{\mathrm{Cov}\left( b,~\frac{a}{b} \right)}{\mathbb{E}[b]}$, we have
Hence we know that $\widetilde{wf(x)}$ is not the unbiased estimator of $\mu_q$. Although in some materials, we may derive that if $N \rightarrow \infty$, $\lim\limits_{N \rightarrow \infty} \widetilde{wf(x)} = \mu_q$, this conclusion could not ensure an unbiased estimation. In some materials, such kind of estimator is described as “asymptotically unbiased”.
To be frank, almost all materials about self-normalized importance sampling do not mention the estimator for variance ($\widetilde{s^2}$), thus we do not talk about it here. First, for any $x_i$, we have
Consider the definition of the standard error, we would have
In this derivation, we make 2 approximations. In the first approximation, we apply $\mathbb{E} \left[\frac{a}{b}\right] = \frac{ \mathbb{E}[a] }{ \mathbb{E}[b] }$, which causes the error brought from the covariance removed. In the second approximation, we know that $\mathbb{E}_p [w] = 1$ while $\frac{1}{N} \mathbb{D}_p [ w ] < \frac{1}{N}$, which means we remove the variance in the dominator. Hence we could only get an approximation of the squared standard error $\widetilde{\varepsilon^2}$ here. If we apply another approximation that $\widetilde{wf(x)} \approx \mu_q$, we could get an estimator
Apparently this is not an unbiased estimator for both $\widetilde{\varepsilon^2}$ and $\varepsilon^2$. If we substitute that $w_i = 1$, we would find that $\widetilde{\mathrm{SE}^2} = \frac{ \sum_{i=1}^N \left(f(x_i) - \overline{f(x)} \right)^2 }{ N^2 }$, which is different from the unbiased estimator that we have derived in $\eqref{fml:mc:unb-se}$.
Simulation results
To show that why the unbiased estimator is so important and how importance sampling works, we write codes to perform Monte-Carlo simulation.
Let $p(\cdot) \sim \mathcal{N}(\mu_n,~\sigma_n)$, $q(\cdot) \sim \chi^2 (k)$. Then we could use different parameters to test different estimators.
For expectation, we test the estimators for direct sampling, importance sampling (where we use $\frac{q(x)}{p(x)}$) and self-normalized importance sampling.
For variance, we compare the estimators for direct sampling and importance sampling. Both estimators are unbiased, however, they are for different variance value.
For standard error, we compare the estimators for 3 methods. Note that the last estimator, i.e. the estimator in self-normalized importance sampling is biased.
For all experiments, we set $k=10$. Then here are the results:
| $(\mu_n,~\sigma_n)$ | Expectation | Variance | Standard error | |
|---|---|---|---|---|
| (20, 10) | ||||
| (9, 6) | ||||
| (20, 3) |
We use dashed red line to show the true value calculated by integration. The blue curves, orange curves and green curves are Monte-Carlo simulation results by direct sampling, by importance sampling and self-normalized importance sampling respectively. Generally we could find that the standard error estimator in self-normalized importance sampling converges to 0 more slowly, because it is not only a biased estimator, but also an estimator for an approximate standard error.
The blue density function is $\chi^2(10)$. Then we know that
In the first experiment, we use a widely spread normal distribution which is not close to the target distribution to perform importance sampling. We could find that the convergence of the expectation is fairly good, while it has an obviously larger variance and standard error compared to direct sampling.
In the second experiment, we adjust the parameters of the normal distribution so that it is closer to the target distribution. The results are even better than the first experiment. Consider that the maximum of the probability density function of the target distribution is about 1.5 times of that of the normal distribution. And we find that the variance of the importance sampling is almost also 1.5 times than the original variance, which proves our theory.
In the last experiment, we use a bad distribution to sampling the target distribution. The distribution is bad because there is a small probability when we get a sample that appears in the target distribution with a high possibility. There is no surprise that the simulation of the importance sampling is in very low quality with an extremely high variance and nearly unconverged standard error.
Metropolis-Hastings algorithm
Check here to see Inverse transform sampling in Wikipedia.
In the above parts, we are keeping discussing about the single value function with a single value distribution. In somehow, some conclusions are also satisfied in the case of multi-variate distribution. However, an inevitable challenge is how could we generate samples from an arbitrary multi-variate distribution.
Consider the case that we only have single value. In this case, we usually use “inverse transform sampling”. The algorithm could be described as
- For the arbitrary probability density function $p(\cdot)$, apply integration and calculate its corresponding cumulative density function $F(\cdot)$.
- Find the inverse function $F^{-1}(\cdot)$.
- Generate a random variable $r \sim U(0,1)$, apply that $x = F^{-1}(r)$, where $x$ is the sample in the distribution $x \sim p(x)$.
Such a technique is not applicable in many cases. For example, if we could not find a inverse of $F(\cdot)$, it may take extra pre-processing to fix this problem. Furthermore, in most cases, we need to deal with multi-variate distribution which makes it very difficult to find a multi-output function $\mathbf{x} = F^{-1}( r )$. Hence we need to find some alternatives to allow us simulate multi-variate sampling. Metropolis-Hastings algorithm is such an approach. Since it is based on Monte-Carlo Markov chain (MCMC), to learn this method, we need to recall the some key points of the theory of Markov chain.
Stationarity of a Markov chain
Consider a Markov chain model. In the example shown in the following figure, we have 5 available states, where the probability that current state is the ith state is $p_i$ (hence the distribution could be represented as $\mathbf{p}$). Denote the probability transition matrix as $\mathbf{Q}$, then we know that the probability that state i turning into state j is $q_{i \rightarrow j}$. We call “a step” in Markov chain
If we have $\mathbf{p}(t+1) = \mathbf{p}(t)$, we call such a distribution as stationary distribution, or $\boldsymbol{\pi}$. Because $\boldsymbol{\pi} = \boldsymbol{\pi} \mathbf{Q}$, for the probability of the jth state, we have
Here we introduce the reversibility condition. If we have such an equation for any i and j,
then the markov chain has the stationary distribution. To understand this condition, we could derive that
Hence we know such a Markov chain is stable. The reversibility condition is a sufficient condition for stationarity. However, a stationary distribution may not require reversibility condition.
Note that if $\mathbf{Q}$ is confirmed, we would have the unique solution for $\boldsymbol{\pi}$. So if we have already known $\mathbf{Q}$ and $\boldsymbol{\pi}$ which could ensure the reversibility condition, then we could confirm that $\boldsymbol{\pi}$ is the steady distribution for this Markov chain.
Continuous form of reversibility condition
Consider that we have a distribution $p(x)$ in continuous form. For any $x,~y$, if we have the transition probability $q(x,~y) := P(x \rightarrow y | x)$ which is continuous, then the reversibility condition could be rewritten as
where $\int_{y \in \Omega} q(x,~y) \mathrm{d} y= 1$. Hence we know that
which ensures the stationarity of the continuous Markov chain.
What is Metropolis-Hastings algorithm
Check here to see Metropolis–Hastings algorithm in Wikipedia.
Metropolis-Hastings algorithm is a MCMC sampling method. The idea of this method is constructing a reversible Markov chain. Suppose that we have a target density function $p(\boldsymbol{\theta})$ which is multi-variate. If we know any $\boldsymbol{\theta}$, we have a jumping distribution $q(\boldsymbol{\theta},~\boldsymbol{\theta}^{\dagger}) := q(\boldsymbol{\theta}^{\dagger} | \boldsymbol{\theta})$. Then given a initial variable $\boldsymbol{\theta}_0$, the algorithm could be described as follows
- For any $t$ ($t \geqslant 0$), generate a sample $\boldsymbol{\theta}^{\dagger} \sim q(\boldsymbol{\theta}_t,~\boldsymbol{\theta}^{\dagger})$.
- Calculate the acceptance rate $\alpha (\boldsymbol{\theta}^{\dagger} | \boldsymbol{\theta}_t) = \max \left( \frac{ p(\boldsymbol{\theta}^{\dagger}) q(\boldsymbol{\theta}^{\dagger},~\boldsymbol{\theta}_t) }{ p(\boldsymbol{\theta}_t) q(\boldsymbol{\theta}_t,~\boldsymbol{\theta}^{\dagger}) },~ 1 \right)$.
- Generate a random variable $r \sim U(0,~1)$. If $r < \alpha$, accept the new sample, i.e. $\boldsymbol{\theta}_{t+1} = \boldsymbol{\theta}^{\dagger}$. Otherwise reject it, i.e. $\boldsymbol{\theta}_{t+1}=\boldsymbol{\theta}_t$.
- Go back to step 1 until we get a group of samples which are unrelated to each other.
In the real application, selecting a good jumping distribution is very important. It is obvious that since $q$ is also a multi-variate distribution, it is necessary to use a simple distribution. Otherwise, we have no approach to generate samples over $q$ quickly. Generally, we use the multi-variate normal distribution, i.e. $q(\cdot|\boldsymbol{\theta}) \sim \mathcal{N}(\boldsymbol{\theta},~\Sigma(\boldsymbol{\theta}))$. We may use different covariance matrices at different position $\boldsymbol{\theta}$. If we make the covariance matrix homogenous, i.e. $\Sigma$ is diagonal, then we have $q(\boldsymbol{\theta}_t,~\boldsymbol{\theta}^{\dagger}) = q(\boldsymbol{\theta}^{\dagger},~\boldsymbol{\theta}_t)$. In this case, the Metropolis-Hastings algorithm degenerates to Metropolis algorithm:
- For any $t$ ($t \geqslant 0$), generate a sample $\boldsymbol{\theta}^{\dagger} \sim q(\boldsymbol{\theta}_t,~\boldsymbol{\theta}^{\dagger})$.
- Calculate the acceptance rate $\alpha (\boldsymbol{\theta}^{\dagger} | \boldsymbol{\theta}_t) = \max \left( \frac{ p(\boldsymbol{\theta}^{\dagger}) }{ p(\boldsymbol{\theta}_t) },~ 1 \right)$.
- Generate a random variable $r \sim U(0,~1)$. If $r < \alpha$, accept the new sample, i.e. $\boldsymbol{\theta}_{t+1} = \boldsymbol{\theta}^{\dagger}$. Otherwise reject it, i.e. $\boldsymbol{\theta}_{t+1}=\boldsymbol{\theta}_t$.
- Go back to step 1 until we get a group of samples which are unrelated to each other.
Stationarity of Metropolis-Hastings algorithm
The continuous transition in Metropolis-Hastings algorithm fulfils the reversibility condition. Note that if a sample is accepted, we have
If the sample is rejected, since the state does not change, we have $\boldsymbol{\theta}^{\dagger} = \boldsymbol{\theta}$, so no matter how much is the probability of the rejection, the reversibility condition still holds.
Estimate the convergence of Metropolis-Hasting algorithm
First, let us introduce the covariance and autocorrelation. For two random variable, the covariance is
Then if we have a sequence $\{ \theta_i \}$ could define the autocorrelation:
Here we introduce an estimator for autocorrelation:
where $\overline{\theta} = \frac{1}{N} \sum_{i=1}^{N} \theta_{i}$.
Actually, $\hat{\rho}_{k}$ is a biased estimator, unless we use the true values of expectation and variance instead of unbiased estimation in $\eqref{fml:mh:rhoest}$. In the following part, we would assume that we could use this estimator to calculate autocorrelation.
AR1 process
To learn how Metropolis-Hastings algorithm converges, we need to generalize this algorithm as a first-order autoregressive process (AR1). This process could be formulated as the following sequence:
where $\varepsilon \sim \mathcal{N} (0,~\sigma)$ is a random noise, $\alpha \in (0,~1)$. Hence we have
where $\varepsilon_1,~ \varepsilon_2,~ \cdots$ are a series random variables which are unrelated to each other.
Hence we could find that $\mathbb{E} \left[ \theta_t - \mu \right] = 0$, i.e.
We could calculate the variance
and the covariance
Hence the truth of the autocorrelation is $\rho_k = \alpha^k$. Since we are using expectation operator to draw this conclusion, $\alpha^k$ is the true value, not an estimator.
Standard error of AR1 process
To learn the standard error, first, we need to calculate
Then we calculate
Consider the covariance matrix
The standard error is the average of all elements in the covariance matrix. Hence we have
Hence we know that $\mathrm{SE} \approx \frac{\sigma}{\sqrt{N}} \sqrt{\frac{1 + \alpha}{1 - \alpha}}$. This is not an estimator, but an approximation to the true value. Since $\rho_1 = \alpha$, we could find that the less autocorrelation in the sequence, the smaller the standard error is.
Simulation results
To test the performance of Metropolis-Hastings algorithm, we design such a probability density function:
| Our specially designed pdf |
|---|
In somehow, this function is not very suitable for applying Metropolis-Hastings, because in most of the input domain, $p(\cdot)=0$. But we could still test it. In this example, we generate a distribution with the mean value of [0.5, 1.0].
We use $\mathbf{x} = \mathbf{x} + \mathcal{N}(0, 0.005)$ to update the sample. Beginning at [0.4, 0.8], after 20000 steps, we choose the last 10000 samples. The results are shown as follows
| Test1 | Test2 | Test3 |
|---|---|---|
We use blue points to show the samples lying in the right domain. The red points represent samples with a $p(\cdot)=0$. The performance of the algorithm is not very stable.
Here we show the estimation of the standard error. Because we have already known the true expectation ([0.5, 1.0]), our estimation for standard error is unbiased here. The result is shown in the following figure:
| Estimation for standard error |
|---|
To prevent the cases that the samples are not in valid domain, we could adapt Metropolis-Hastings algorithm as below
- For any $t$ ($t \geqslant 0$), generate a sample $\boldsymbol{\theta}^{\dagger} \sim q(\boldsymbol{\theta}_t,~\boldsymbol{\theta}^{\dagger})$.
- Calculate the acceptance rate $\alpha (\boldsymbol{\theta}^{\dagger} | \boldsymbol{\theta}_t) = \max \left( \frac{ p(\boldsymbol{\theta}^{\dagger}) q(\boldsymbol{\theta}^{\dagger},~\boldsymbol{\theta}_t) }{ p(\boldsymbol{\theta}_t) q(\boldsymbol{\theta}_t,~\boldsymbol{\theta}^{\dagger}) },~ 1 \right)$.
- If $p(\boldsymbol{\theta}^{\dagger}) < T$, where $T$ is a threshold, let $\alpha=0$
- Generate a random variable $r \sim U(0,~1)$. If $r < \alpha$, accept the new sample, i.e. $\boldsymbol{\theta}_{t+1} = \boldsymbol{\theta}^{\dagger}$. Otherwise reject it, i.e. $\boldsymbol{\theta}_{t+1}=\boldsymbol{\theta}_t$.
- Go back to step 1 until we get a group of samples which are unrelated to each other.
By improving the algorithm, we could set a larger $\sigma$, i.e. $\mathbf{x} = \mathbf{x} + \mathcal{N}(0, 1)$. Then we would have new results:
| Sampling iterations | Test Results (Last 10K samples) |
|---|---|
| 20K | |
| 100K | |
| 1M |
| Sampling iterations | Estimation for standard error (improved algorithm) |
|---|---|
| 20K | |
| 100K | |
| 1M |
We could find that after improving the algorithm, there is no point lying out of the valid domain. Furthermore, it seems that in these simulations, the autocorrelation does not influence the standard error very much, which may be caused by our distribution.
Gibbs sampling
We have known that although Metropolis-Hastings algorithm could be a good simulation, the results may be influenced by the distribution of jumping probability. To get rid of such influence, we may let $q(\cdot) \sim U(\mathbf{x} \in \Omega)$, which means using uniform distribution in the whole input domain. However, if we could generate the conditional sampling of the target distribution $p(\cdot)$, Gibbs sampling would be a better choice.
What is Gibbs sampling
Gibbs sampling could be viewed as a special case of the Metropolis-Hastings algorithm. It is also used to generate samples for a multi-variate distribution. Consider that we have a multi-variate variable $\boldsymbol{\theta} \sim p(\boldsymbol{\theta})$, If we set the initial guess as $\boldsymbol{\theta}^{(0)}$, then Gibbs sampling could be described as
For any $k \geqslant 0$, the kth iteration of Gibbs sampling could be described as
- $\theta^{(k+1)}_1 \sim p\left(\theta_1 \left| \theta_2=\theta^{(k)}_2,~ \theta_3=\theta^{(k)}_3,~ \theta_n=\theta^{(k)}_n \right.\right)$
- $\theta^{(k+1)}_2 \sim p\left(\theta_2 \left| \theta_1=\theta^{(k+1)}_1,~ \theta_3=\theta^{(k)}_3,~ \theta_n=\theta^{(k)}_n \right.\right)$
- $\theta^{(k+1)}_3 \sim p\left(\theta_3 \left| \theta_1=\theta^{(k+1)}_2,~ \theta_2=\theta^{(k+1)}_2,~ \theta_n=\theta^{(k)}_n \right.\right)$
- ...
- $\theta^{(k+1)}_n \sim p\left(\theta_n \left| \theta_1=\theta^{(k+1)}_1,~ \theta_2=\theta^{(k+1)}_2,~ \theta_n=\theta^{(k+1)}_{n-1} \right.\right)$
Although it is much easier for us to sample from a conditional distribution with a single variable and other parameters fixed, it may still take us some efforts to calculate the conditional sampling. So Gibbs sampling may be not efficient in some cases that even the conditional sampling is difficult.
Staionarity of Gibbs sampling
Consider that we are running the ith step of Gibbs sampling. In this case, we know that
The acceptance rate is
which means as Gibbs sampling always accept the new sample.
Because we could find that in any step of any iteration, $q(\boldsymbol{\theta},~\boldsymbol{\theta}^{\dagger}) = p(\boldsymbol{\theta}^{\dagger})$, and $q(\boldsymbol{\theta}^{\dagger},~\boldsymbol{\theta}) = p(\boldsymbol{\theta})$, the reversible condition is satisfied obviously.
Approximate the marginal density function
If we could get a sampling result, then we could approximate the marginal distribution by Monte-Carlo method. In each iteration, we could calculate that
We could select $M$ samples from the Gibbs sampling results randomly. For each sample, calculate the continuous probability density function for $\theta_i$, then use Monte-Carlo method to calculate the expectation of all other variables. Then the result is the marginal distribution.
Estimate the standard error of a sampling method
Denote that $\mathrm{Cov}(\theta_i,~\theta_{i+k}) = \gamma(k)$. From the previous derivation for AR1 process, we could get the similar conclusion that
Consider the covariance matrix
Hence
Simulation results
To test the performance of Gibbs sampling, we still use the same density function which is used in testing Metropolis-Hastings algorithm:
| Our specially designed pdf |
|---|
To perform Gibbs sampling, we need to specify the conditional probability. Denote the center of the distribution as $[\mu_x,~\mu_y]$, and the stretch of the distribution as $[S_x,~S_y]$. The parameters are shown in the following figure.
| A schema for our 2D distribution |
|---|
For a known $x_0$ or $y_0$, we have
where $\mathrm{Triangle}( x | \mu,~S )$ is the 1D triangular distribution. The center locates in $\mu$ while the boarder which defines $p(x) = 0$ is $[ \mu-S,~\mu+S ]$.
We still begin with [0.4, 0.8]. After 20000 steps, we choose the last 10000 samples. The results are shown as follows
| Test results |
|---|
It seems that the performance is very well.
Here we show the estimation of the standard error. Because we have already known the true expectation ([0.5, 1.0]), our estimation for standard error is unbiased here. The result is shown in the following figure:
| Estimation for standard error |
|---|
We could use the unbiased estimator $\overline{\mathbf{x}}$ to estimate the expectation, then we could define the distance from the estimator to the real value, i.e. $\lVert \overline{\mathbf{x}} - \boldsymbol{\mu} \rVert_2$. The results are shown as follows.
| Estimation for distance from average to true expectation |
|---|
Since we do not change the estimation size $N$ but only change the start step of the estimation. Theoretically both the estimation for standard error and the estimation for distance error should be stable. The results prove the theory.
Source code
The source code of this inspection has been uploaded to Github, to learn more about the codes, please check the following link:


She sights a Bird—she chuckles—
She flattens—then she crawls—
She runs without the look of feet—
Her eyes increase to Balls—
Her Jaws stir—twitching—hungry—
Her Teeth can hardly stand—
She leaps, but Robin leaped the first—
Ah, Pussy, of the Sand,
The Hopes so juicy ripening—
You almost bathed your Tongue—
When Bliss disclosed a hundred Toes—
And fled with every one—


Dashboard
Summary
Total:
Running:
Idle:
Operation Rate %:
Status
| # | APP Name | Time |
|---|
Errors
| # | APP Name | Date | Time | Error |
|---|
$app_name$ $app_status$
System Memory (MB)
Total:
Free:
Used:
Used %:
APP Service Memory (MB)
Allocated:
Free:
Used:
Used %:
Test
| # | Date | Time | Result |
|---|
Error
| # | Date | Time | Error |
|---|
API Time
| # | Path | Time (ms) |
|---|


--- layout: title: Use cases and productions locale: en permalink: "/en/blog" --- {% include head.html %}
{{ref.article-description}}


Cicada is distributed as either a set of C or C++ source files. The C version of Cicada is designed to embed (or be embedded in) the user’s C code, whereas the C++ version merges with a user’s C++ code. See the previous section for details on how to do do this.
The simplest way to build the executable is to use the command-line tool make along with a C/C++ compiler such as gcc; free versions of these exist for many platforms. To build Cicada using make, put everything---source files, header files, Makefile and ‘.cicada’ files---into the same directory; navigate to that directory from a command prompt; and then type “make cicada CC=gcc” for a C program, or just “make cicada” for a C++ program. This should produce a cicada executable in that same directory.
Once Cicada has successfully been compiled, we can run it by opening up a command terminal, navigating into Cicada’s directory and typing “./cicada” (UNIX and Mac), or just “cicada” (PC/Windows). Cicada should present a command prompt, looking like this:
>
Last update: July 25, 2021


Highly On: the Wikipedia article would have contributed to constructing knowledge of the webpage./ I would read this Wikipedia article to learn more about the webpage.
Slightly On: this Wikipedia article covers a certain topic of the webpage, but it is not the main topic.
Slightly Off: this Wikipedia article is somewhat related to the webpage in a broad sense, but not directly related.
Highly Off: this Wikipedia article has nothing to do with the webpage.


5 August 2005
The Ministry of Health in Viet Nam has today confirmed an additional three cases of human infection with H5N1 avian influenza. One case was reported in Ha Tay Province, one in Tra Vinh Province, and one in Ho Chi Minh City. The patients from Tra Vinh and Ho Chi Minh City died.
The newly confirmed cases bring the total in Viet Nam since mid-December 2004 to 63 cases, of which 20 were fatal.


{% extends "base_admin.html" %} {% block content_admin %}
This page is for the management of App Flask & Admin


"""
Contains constants relating to robot configuration; for example, Talon CAN IDs
and frame dimensions.
"""
import wpilib
# Teleop control constants. Can be loaded from Preferences.
fwdAxis = 1 #: Forward/Backward axis
strAxis = 0 #: Left/Right axis
rcwAxis = 2 #: Rotation axis
liftAxis = 2 #: Lift control axis on throttle
fwdInv = True #: Fwd/Bwd axis inverted
strInv = True #: L/R axis inverted
rcwInv = True #: Rot axis inverted
teleop_speed = 370
turn_sensitivity = 0.25
lift_deadband = 0.25 # deadband
# Wraps the Preferences API to provide an alternative to all of the
# getInt/getString/getWhatever methods
def __load_preference(key, backup):
prefs = wpilib.Preferences.getInstance()
getMethod = None
putMethod = None
if isinstance(backup, str):
getMethod = prefs.getString
putMethod = prefs.putString
elif isinstance(backup, bool):
getMethod = prefs.getBoolean
putMethod = prefs.putBoolean
elif isinstance(backup, int):
getMethod = lambda k, b: int(prefs.getInt(k, b)) # noqa: E731
putMethod = lambda k, v: prefs.putInt(k, int(v)) # noqa: E731
elif isinstance(backup, float):
getMethod = lambda k, b: float(prefs.getFloat(k, b)) # noqa: E731
putMethod = lambda k, v: prefs.putFloat(k, float(v)) # noqa: E731
if not prefs.containsKey(key):
putMethod(key, backup)
return backup
else:
return getMethod(key, backup)
[docs]def load_control_config():
"""
Load configurable constants using the Robot Preferences API.
Do not call this at module level (otherwise it might try to access parts of
WPILib before they have been initialized).
"""
global fwdAxis, fwdInv, strAxis, strInv, rcwAxis, rcwInv, teleop_speed
global turn_sensitivity, liftAxis, liftInv, lift_deadband
fwdAxis = __load_preference('Control: Forward-Backward Axis', backup=1)
fwdInv = __load_preference('Control: Fwd-Bwd Axis Inverted', backup=True)
strAxis = __load_preference('Control: Left-Right Axis', backup=0)
strInv = __load_preference('Control: L-R Axis Inverted', backup=True)
rcwAxis = __load_preference('Control: Rotation Axis', backup=2)
rcwInv = __load_preference('Control: Rot Axis Inverted', backup=True)
liftAxis = __load_preference('Control: Lift Control Axis', backup=3)
liftInv = __load_preference('Control: Lift Control Inverted', backup=True)
lift_deadband = __load_preference(
'Control: Lift Control Deadband', backup=0.25
)
teleop_speed = __load_preference('Control: Teleop Speed', backup=370)
turn_sensitivity = __load_preference(
'Control: Turn Sensitivity', backup=0.25
)
#: Swerve module hardware configuration.
#: List of tuples of form ('module name', steer_id, drive_id)
#: See swerve/swerve_drive.py
swerve_config = [
('Back Right', 6, 4),
('Back Left', 7, 5),
('Front Right', 8, 9),
('Front Left', 11, 10),
]
#: Lift motor contorller CAN IDs. Currently dummy values.
lift_ids = {
'left': 31,
'right': 42
}
# Claw motor controller CAN ID(s).
claw_id = 2
# The length of the chassis (units do not matter as long as they match)
chassis_length = 24
# The width of the chassis (units do not matter as long as they match)
chassis_width = 27
# Winch Motor CAN ID
winch_id = 35


×
Reactive CEP Example
Distance Pattern Type
Min Distance Pattern
Max Distance Pattern
Average Distance Pattern
Window
Tumbling Count Window
Sliding Count Window
Hopping Count Window
Tumbling Temporal Window
Hopping Temporal Window
Location
help_outline
Latitude:
Longitude:
Pattern Assertion
New Event Type Id
Set Pattern
Unset Pattern
View Code
Event Type Id
Detection Time
JSON


×
Notes
Reminder
Labels
Create New Label
Archive
Trash
Labels
Settings
Send Feedback
Help
App Download
Keyboard Shortcuts
Privacy
Terms
Open Source Licenses
☰
Google Keep


## Mako can find this base file because it's in the settings.py -> MAKO_TEMPLATES_DIRS list <%inherit file="base.htm" /> <%block name="title">


{{job.name}} # Invoice Item Name Price #{{expense.ticketNum}} {{expense.item}} {{expense.price | currency}} Submit


| Reference | Class List | Member Functions | Typedefs | File List |
This container implements a scheduling policy based on task priorities. The task with highest priority will be the first to be removed. It must be possible to compare two tasks using operator<.
| Task | A function object which implements the operator() and operator<. operator< must be a partial ordering. |
Definition at line 193 of file scheduling_policies.hpp.
Public Types | |
| typedef Task | task_type |
Public Member Functions | |
| void | clear () |
| bool | empty () const |
| void | pop () |
| bool | push (task_type const &task) |
| size_t | size () const |
| task_type const & | top () const |
Protected Attributes | |
| std::priority_queue< task_type > | m_container |
|
|
Indicates the scheduler's task type.
Definition at line 196 of file scheduling_policies.hpp. |
|
|
Removes all tasks from the scheduler. Definition at line 248 of file scheduling_policies.hpp. |
|
|
Checks if the scheduler is empty.
Definition at line 241 of file scheduling_policies.hpp. |
|
|
Removes the task which should be executed next. Definition at line 215 of file scheduling_policies.hpp. |
|
|
Adds a new task to the scheduler.
Definition at line 207 of file scheduling_policies.hpp. |
|
|
Gets the current number of tasks in the scheduler.
Definition at line 232 of file scheduling_policies.hpp. |
|
|
Gets the task which should be executed next.
Definition at line 223 of file scheduling_policies.hpp. |
|
|
Internal task container.
Definition at line 199 of file scheduling_policies.hpp. |


Basic fs-person-gender demo
Editable fs-person-gender demo
Editing fs-person-gender demo
Basic fs-person-gender-selector demo
Setting the gender on fs-person-gender-selector demo


Ant McPartlin's ex-wife Lisa Armstrong makes new dig with cosmetic range name
The make-up artist has been locked in a divorce battle with the TV presenter.
by Dailyrecord.co.ukThe ex-wife of Ant McPartlin has taken a swipe at the TV host with the name of her new cosmetic range.
Lisa Armstrong has reportedly called the products after their beloved pet, who they are locked in a custody batter over.
According the The Sun, the make-up artist will name her range Lisa Love Hurley.
The 43-year-old has been in a divorce battle with Ant as she wants to keep dog Hurley full-time - but Ant is determined to spend half the time with the pup.
“It’s a way for her to say ‘F*** you’ to Ant, and a declaration of how much she loves Hurley," a source told the Mirror .
"He’s been there by her side through all the mess Ant created.
“She has always been upset that Ant couldn’t just let her have custody of Hurley who she sees as her baby.
“She knows that Ant loves Hurley as well but she really hoped he would be the bigger man and let her at least have custody of him after all she has been through."
It comes as Lisa has rejected a proposed settlement of £31million - half of Ant's £62million fortune - as well as their former marital home in London's upmarket Chiswick.
Instead, Lisa's said to be wanting her day in court and has so far refused to sign Ant's non-disclosure agreement (NDA) that would be part of any financial settlement her legal team thrashes out with his.


Arlo Pro 3 review: Two steps forward, one step back
by Stephen SchenckWireless cameras have had a huge impact on how we're able to keep a watchful eye over our homes, possessions, and loved ones. Modern systems are as feature packed as they are easy to install, and the field is still new enough that we're seeing real improvements generation-to-generation. Earlier this year, Arlo launched a 4K option in the form of its Ultra cam, and now we're getting a fresh entry in its Pro series with the release of the Arlo Pro 3.
The Pro 2 was already a fantastic wireless camera, with a solid 1080p resolution, battery life extending into weeks, and free 7-day cloud storage for recorded clips. Now the Pro 3's out to try and one-up its predecessor, offering a new 2K resolution, integrated spotlight, wide-angle recording, and more. But for all those changes, does the Pro 3 ultimately offer a better user experience? Let's check it out.
Specs
| Camera | 2560 x 1440 2K-class |
| Field of view | 160 degrees |
| Night vision | IR LEDs for B&W, white LED "spotlight" for color |
| Audio | Microphone array with noise and echo cancellation |
| Power | 4,800 mAh battery |
| Connectivity | Dedicated hub, 802.11 b/g/n |
| Local storage | USB support for FAT32 drive (up to 2TB) |
| Cloud storage | None free, available as paid add-on ($3/mo./camera, $10 for up to 10) |
The Good
| Video quality | A distinct improvement over the 1080p Pro 2: crisper details, wider angles. |
| Rugged design | Feels a lot more resilient and weatherproof than last-gen. |
| Installation | From hanging camera mounts to configuring software, Arlo makes getting set up a breeze. |
The Not So Good
| Price | We knew it was going to be expensive, but $500 for a two-pack is still a lot to stomach. |
| Arlo Smart | Delivering object recognition as software-as-a-service and killing off free 7-day cloud recording makes a Smart subscription practically mandatory. |
| Charging | Arlo was so close to nailing charging on the Pro 3, but a proprietary cable is a step in the wrong direction. |
Design, hardware, what's in the box
The Arlo Pro 3 goes with the same basic design as the Ultra, and compared to the Pro 2, this is a big improvement. Rather than giving us another camera body peppered with hinges, buttons, and flaps, the Pro 3 encloses its hardware in a sleek, well-protected shell. The camera body pops out of the frame to access the battery, but is normally held snugly in place with a gasket to keep all the forces of nature at bay.
Arlo includes a pair of mounting options for the camera two-pack I checked out. The standard magnetic mount gets an upgrade this year, moving from a convex to concave design. This gives you a bit more freedom in setting up your viewing angle, but mostly just helps the cameras themselves feel more durable. We also get a pole mount, which works great when installing the Pro 3 on a ceiling.
While the pole mount is arguably more secure, the magnetic mount easily installs with just a single screw.
My only real complaint with the new look is that charging is now done with a proprietary magnetic cable. I certainly don't miss the rubber flap covering the Pro 2's micro-USB port, but while more robust, this seems destined to get lost (it's $25 to replace). The new camera-shell design also seems like it would be ideal to let you pop the camera out of its frame to charge, while not forcing you to reconfigure viewing angles later. Unfortunately, you still have to pull the whole thing down, as the magnetic adapter only works with the shell present.
Just make sure you keep the proprietary magnetic charging adapter in a safe place.
The big news for the Pro 3's camera hardware is the move to 2K resolution: 2560 x 1440. We also see viewing angles expand from 130 to 160 degrees, and the previous B&W IR-only night-vision mode joined by an option for full color thanks to the presence of an LED spotlight. We'll talk about image quality and operation in just a moment.
Arlo's base station (formally a "SmartHub") is just slightly different this year, dropping one of its USB ports for external storage (it's still got the one), and offering fewer status LEDs. Probably the most significant move there is that the system's siren relocates from within the base station itself to being a part of each individual Pro 3 camera.
A more streamlined hub still offers USB support for external storage, but drops its siren.
Video and app features
Well, the move to 2K has really paid off: Video shot with the Pro 3 looks great.  Wide-angle filming, while not quite as wide as on the Ultra, is a great asset, and you can always tighten the field of view up if you don't need to grab all the corners. There's even an optional auto-tracking mode that digitally zooms in to follow the action, though I prefer to keep the whole scene in view.
Motion detection works well, and is highly configurable with sensitivity options and customizable zones. If you actually want video clips of all the chipmunks running around your yard (like me), the Pro 3 is happy to comply, or can easily be set to only look for more substantial motion.
Speaking of that, advanced features like person, animal, vehicle, and package detection are available, helping to curate the motion events that need your attention. Unfortunately, all this is only yours as part of the Arlo Smart subscription service, though you do get a three-month free trial.
Camera overview, Smart Alert settings, videos classified by type (left to right)
That's also the only way to get cloud storage for the Pro 3, which drops the free 7-day rolling storage from the Pro 2. Honestly, this is an absolute deal-breaker for me personally, and with how expensive Arlo hardware is, this just feels like nickel-and-diming. We're not talking a lot of storage, either: The average for my motion-triggered clips over a 7-day period comes to under 1GB per camera. A $200 camera should be able to come with a measly gigabyte of cloud storage.
This might be less of a problem if the Arlo hub's external USB storage were worth a damn, but there's no way to watch those recordings from the app itself; you have to physically move your USB device over to a computer to view. Even the bargain-basement $20 Wyze Cam lets you stream its local microSD storage through the app. This feels less like a technical limitation and more a move to generate subscribers for Arlo Smart. Either way, it's half baked and makes me cross.
The Arlo app itself hasn't changed a ton since I got to know it with the Pro 2, with a healthy assortment of options to dive into, and basic rule support — motion on one camera can trigger another to start recording. While familiar and functional, some of the same problems I've seen before remain, like lengthy delays when manually connecting to cameras for live views, and the occasional duplicate recording.
Should you buy it?
Arlo Pro 3
8/10
Maybe. The Pro 3 camera hardware is an undeniable improvement over the Pro 2, and while I feel Arlo didn't go far enough in a few areas, and missed some opportunities to make things even better, I'm very satisfied with the overall direction. But while the Pro 2 was a great investment, to get any sort of reasonable operation out of the Pro 3 you'll need to pay between $70 and $120 a year for Arlo Smart and its cloud connectivity. Arlo's always catered to shoppers willing to spend top-dollar, and the Pro 3 only further cements that role.
If price is no object, you should be very happy with the Arlo Pro 3. The two-pack system reviewed here lists for $500, with add-on cameras at $200 a pop (though it's easy to save with bundle deals). But if you'd rather not worry about ongoing service fees, there are plenty of other options out there, and it's as good a time as any to pick up the Pro 2 for cheap.


2 Map Neighborhoods
When considering the health of persons, we have to also consider the neighborhood environment. Sometimes this is looking at neighborhood level health outcomes, like premature mortality at the census tract scale, or cumulative COVID rates by zip code. Sometimes we’re interested in neighborhood factors like poverty, access to affordable housing, or distance to nearest health provider, or pollution-emitting facility. These measurements of the “social determinants of health” at the neighborhood scale are increasingly urgent in modern public health thinking, and are thought to drive and/or reinforce racial, social, and spatial inequity
In this module, we’ll learn about the basics of thematic mapping – known as choropleth mapping – to visualize neighborhood level health phenomena. This will allow you to begin the process of exploratory spatial data analysis and hypothesis generation & refinement.
2.1 Clean Attribute Data
Let’s consider COVID-19 cases by zip code in Chicago. We’ll upload and inspect a summary of cases from the Chicago Data Portal first:
## ZIP.Code Week.Number Week.Start Week.End Cases...Weekly
## 1 60603 39 09/20/2020 09/26/2020 0
## 2 60604 39 09/20/2020 09/26/2020 0
## 3 60611 16 04/12/2020 04/18/2020 8
## 4 60611 15 04/05/2020 04/11/2020 7
## 5 60615 11 03/08/2020 03/14/2020 NA
## 6 60603 10 03/01/2020 03/07/2020 NA
## Cases...Cumulative Case.Rate...Weekly Case.Rate...Cumulative
## 1 13 0 1107.3
## 2 31 0 3964.2
## 3 72 25 222.0
## 4 64 22 197.4
## 5 NA NA NA
## 6 NA NA NA
## Tests...Weekly Tests...Cumulative Test.Rate...Weekly
## 1 25 327 2130
## 2 12 339 1534
## 3 101 450 312
## 4 59 349 182
## 5 6 9 14
## 6 0 0 0
## Test.Rate...Cumulative Percent.Tested.Positive...Weekly
## 1 27853.5 0.0
## 2 43350.4 0.0
## 3 1387.8 0.1
## 4 1076.3 0.1
## 5 21.7 NA
## 6 0.0 NA
## Percent.Tested.Positive...Cumulative Deaths...Weekly
## 1 0.0 0
## 2 0.1 0
## 3 0.2 0
## 4 0.2 0
## 5 NA 0
## 6 NA 0
## Deaths...Cumulative Death.Rate...Weekly Death.Rate...Cumulative
## 1 0 0 0
## 2 0 0 0
## 3 0 0 0
## 4 0 0 0
## 5 0 0 0
## 6 0 0 0
## Population Row.ID ZIP.Code.Location
## 1 1174 60603-39 POINT (-87.625473 41.880112)
## 2 782 60604-39 POINT (-87.629029 41.878153)
## 3 32426 60611-16 POINT (-87.620291 41.894734)
## 4 32426 60611-15 POINT (-87.620291 41.894734)
## 5 41563 60615-11 POINT (-87.602725 41.801993)
## 6 1174 60603-10 POINT (-87.625473 41.880112)
Each row corresponds to a zip code at a different week. This data thus exists as a “long” format, which doesn’t work for spatial analysis. We need to convert to “wide” format, or at the very least, ensure that each zip code corresponds to one row.
To simplify, let’s identify the last week of the dataset, and then subset the data frame to only show that week. We will be interested in the cumulative case rate. Following is one way of doing this – can you think of another way? Try out different approaches of reshaping data to test your R and “tidy” skills.
## [1] 1 31
## ZIP.Code Week.Number Week.Start Week.End Cases...Weekly
## 1 60603 39 09/20/2020 09/26/2020 0
## 2 60604 39 09/20/2020 09/26/2020 0
## 36 60601 39 09/20/2020 09/26/2020 8
## 37 60602 39 09/20/2020 09/26/2020 0
## 41 60605 39 09/20/2020 09/26/2020 12
## 66 60610 39 09/20/2020 09/26/2020 35
## Cases...Cumulative Case.Rate...Weekly Case.Rate...Cumulative
## 1 13 0 1107.3
## 2 31 0 3964.2
## 36 213 54 1451.4
## 37 21 0 1688.1
## 41 391 44 1420.8
## 66 666 90 1706.9
## Tests...Weekly Tests...Cumulative Test.Rate...Weekly
## 1 25 327 2130
## 2 12 339 1534
## 36 202 4304 1376
## 37 27 460 2170
## 41 291 7160 1058
## 66 500 10680 1281
## Test.Rate...Cumulative Percent.Tested.Positive...Weekly
## 1 27853.5 0.0
## 2 43350.4 0.0
## 36 29328.8 0.0
## 37 36977.5 0.0
## 41 26018.4 0.0
## 66 27371.3 0.1
## Percent.Tested.Positive...Cumulative Deaths...Weekly
## 1 0.0 0
## 2 0.1 0
## 36 0.0 1
## 37 0.0 0
## 41 0.1 1
## 66 0.1 0
## Deaths...Cumulative Death.Rate...Weekly Death.Rate...Cumulative
## 1 0 0.0 0.0
## 2 0 0.0 0.0
## 36 6 6.8 40.9
## 37 0 0.0 0.0
## 41 3 3.6 10.9
## 66 10 0.0 25.6
## Population Row.ID ZIP.Code.Location
## 1 1174 60603-39 POINT (-87.625473 41.880112)
## 2 782 60604-39 POINT (-87.629029 41.878153)
## 36 14675 60601-39 POINT (-87.622844 41.886262)
## 37 1244 60602-39 POINT (-87.628309 41.883136)
## 41 27519 60605-39 POINT (-87.623449 41.867824)
## 66 39019 60610-39 POINT (-87.63581 41.90455)
To clean our data a bit, we’ll just keep the zip code name, and cumulative case rate for the week of September 20th, 2020.
## ZIP.Code Case.Rate...Cumulative
## 1 60603 1107.3
## 2 60604 3964.2
## 36 60601 1451.4
## 37 60602 1688.1
## 41 60605 1420.8
## 66 60610 1706.9
2.2 Merge Spatial Data
Next, let’s merge this data to our zip code master spatial file. Reload if necessary:
## Reading layer `geo_export_54bc15d8-5ef5-40e4-8f72-bb0c6dbac9a5' from data source
## `/Users/maryniakolak/code/Intro2RSpatialMed/data/geo_export_54bc15d8-5ef5-40e4-8f72-bb0c6dbac9a5.shp'
## using driver `ESRI Shapefile'
## Simple feature collection with 61 features and 4 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: -87.94011 ymin: 41.64454 xmax: -87.52414 ymax: 42.02304
## CRS: 4326
## Simple feature collection with 6 features and 4 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: -87.80649 ymin: 41.88747 xmax: -87.59852 ymax: 41.93228
## CRS: 4326
## objectid shape_area shape_len zip geometry
## 1 33 106052287 42720.04 60647 MULTIPOLYGON (((-87.67762 4...
## 2 34 127476051 48103.78 60639 MULTIPOLYGON (((-87.72683 4...
## 3 35 45069038 27288.61 60707 MULTIPOLYGON (((-87.785 41....
## 4 36 70853834 42527.99 60622 MULTIPOLYGON (((-87.66707 4...
## 5 37 99039621 47970.14 60651 MULTIPOLYGON (((-87.70656 4...
## 6 38 23506056 34689.35 60611 MULTIPOLYGON (((-87.61401 4...
Next, merge on zip code ID. The key in the Chi_Zips object is zip
, whereas the key for the COVID data is ZIP.code
. Always merge non-spatial to spatial data, not the other way around. Think of the spatial file as your master file that you will continue to add on to…
#Chi_Zipsf <- merge(Chi_Zips, COVID.39f, by.x = "zip", by.y = "ZIP.Code", all = TRUE)
Chi_Zipsf <- merge(Chi_Zips, COVID.39f, by.x = "zip", by.y = "ZIP.Code")
head(Chi_Zipsf)
## Simple feature collection with 6 features and 5 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: -87.63999 ymin: 41.85317 xmax: -87.60246 ymax: 41.88913
## CRS: 4326
## zip objectid shape_area shape_len Case.Rate...Cumulative
## 1 60601 27 9166246 19804.58 1451.4
## 2 60602 26 4847125 14448.17 1688.1
## 3 60603 19 4560229 13672.68 1107.3
## 4 60604 48 4294902 12245.81 3964.2
## 5 60605 20 36301276 37973.35 1420.8
## 6 60606 31 6766411 12040.44 2289.6
## geometry
## 1 MULTIPOLYGON (((-87.62271 4...
## 2 MULTIPOLYGON (((-87.60997 4...
## 3 MULTIPOLYGON (((-87.61633 4...
## 4 MULTIPOLYGON (((-87.63376 4...
## 5 MULTIPOLYGON (((-87.62064 4...
## 6 MULTIPOLYGON (((-87.63397 4...
2.3 Quantile Maps
Starting with a “classic epi” approach, let’s look at case rates as quantiles. We use the tmap
library, and update the choropleth data classification using the style
parameter. We use the Blue-Purple palette, or BuPu
, from Colorbrewer.
Colorbrewer Tip: To display all Colorbrewer palette options, load the RColorBrewer
library and run display.brewer.all()
– or just Google “R Colorbrewer palettes.”
library(tmap)
tmap_mode("plot")
tm_shape(Chi_Zipsf ) +
tm_polygons("Case.Rate...Cumulative",
style="quantile", pal="BuPu",
title = "COVID Case Rate")
Let’s try tertiles:
2.4 Standard Deviation Maps
While quantiles are a nice start, let’s classify using a standard deviation map. Standard deviation is a statistical technique type of map based on how much the data differs from the mean.
2.5 Jenks Maps
Another approach of data classification is natural breaks, or jenks. This approach looks for “natural breaks” in the data using a univariate clustering algorithm.
tm_shape(Chi_Zipsf ) +
tm_polygons("Case.Rate...Cumulative",
style="jenks", pal="BuPu",
title = "COVID Case Rate")
The first bin doesn’t seem very intuitive. Let’s try 4 bins instead of 5 by changing the n
parameter. In this version, we’ll also had a histogram and scale bar, and move the legend outside the frame to make it easier to view.
2.6 Integrate More Data
To explore potential disparities in COVID health outcomes, let’s bring in pre-cleaned demographic, racial, and ethnic data from the Opioid Environment Policy Scan database. This data is orginally sourced from the American Community Survey 2018 5-year estimate, which you could also pull using the tidycensus
.
## ZCTA year totPopE whiteP blackP amIndP asianP pacIsP otherP hispP
## 1 35004 2018 11762 84.39 13.09 0.00 0.94 0.00 1.57 0.94
## 2 35005 2018 7528 55.22 42.44 0.64 0.00 0.15 1.55 1.37
## 3 35006 2018 2927 96.04 3.21 0.27 0.00 0.00 0.48 0.00
## 4 35007 2018 26328 73.83 13.75 0.04 1.33 0.02 11.01 11.11
## 5 35010 2018 20625 63.07 32.43 0.39 0.65 0.00 3.45 4.10
## 6 35013 2018 40 100.00 0.00 0.00 0.00 0.00 0.00 100.00
## noHSP age0_4 age5_14 age15_19 age20_24 age15_44 age45_49 age50_54
## 1 5.52 787 1950<PHONE_NUMBER> 662 541
## 2 17.48 511 1055<PHONE_NUMBER> 580 469
## 3 14.44 161 413 141 203 878 129 193
## 4 12.41 1891 4161 1619 1400 9947 1993 2067
## 5 22.00 1013 2647 1383 1087 7036 1418 1545
## 6 100.00 0 0 0 0 13 8 19
## age55_59 age60_64 ageOv65 ageOv18 age18_64 a15_24P und45P ovr65P
## 1<PHONE_NUMBER> 8820 7158 10.23 61.97 14.13
## 2<PHONE_NUMBER> 5691 4319 9.72 53.07 18.23
## 3 316<PHONE_NUMBER> 1749 11.75 49.61 19.10
## 4 1713 1315 3241 19178 15937 11.47 60.77 12.31
## 5 1510 1341 4115 16142 12027 11.98 51.86 19.95
## 6 0 0 0 40 40 0.00 32.50 0.00
## disbP
## 1 12.7
## 2 23.2
## 3 20.9
## 4 13.5
## 5 19.6
## 6 0.0
Merge to our master Zip Code dataset.
## Simple feature collection with 6 features and 31 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: -87.63999 ymin: 41.85317 xmax: -87.60246 ymax: 41.88913
## CRS: 4326
## zip objectid shape_area shape_len Case.Rate...Cumulative year
## 1 60601 27 9166246 19804.58 1451.4 2018
## 2 60602 26 4847125 14448.17 1688.1 2018
## 3 60603 19 4560229 13672.68 1107.3 2018
## 4 60604 48 4294902 12245.81 3964.2 2018
## 5 60605 20 36301276 37973.35 1420.8 2018
## 6 60606 31 6766411 12040.44 2289.6 2018
## totPopE whiteP blackP amIndP asianP pacIsP otherP hispP noHSP age0_4
## 1 14675 74.17 5.57 0.45 18.00 0.00 1.81 8.68 0.00 550
## 2 1244 68.17 3.78 5.31 19.45 0.00 3.30 6.51 0.00 61
## 3 1174 63.46 3.24 0.00 27.60 0.00 5.71 9.80 0.00 13
## 4 782 63.43 5.63 0.00 29.67 0.00 1.28 4.35 0.00 12
## 5 27519 61.20 17.18 0.18 16.10 0.03 5.31 5.84 2.39 837
## 6 3101 72.75 2.35 0.00 18.09 0.00 6.80 6.29 0.73 57
## age5_14 age15_19 age20_24 age15_44 age45_49 age50_54 age55_59
## 1 156<PHONE_NUMBER> 976 1009 324
## 2 87 18 91 987 46 53 0
## 3 43 179 172 684 75 47 150
## 4 7 52 168 450 27 47 54
## 5 1279 2172 2282 16364 1766 1520 1824
## 6 44 0 139 1863 213 153 168
## age60_64 ageOv65 ageOv18 age18_64 a15_24P und45P ovr65P disbP
## 1 859 2075 13855 11780 12.37 64.27 14.14 6.4
## 2 5 5 1095 1090 8.76 91.24 0.40 0.2
## 3 50 112 1118 1006 29.90 63.03 9.54 7.3
## 4 92 93 744 651 28.13 59.97 11.89 4.1
## 5 1360 2569 25259 22690 16.19 67.15 9.34 5.3
## 6<PHONE_NUMBER> 2569 4.48 63.33 13.90 1.9
## geometry
## 1 MULTIPOLYGON (((-87.62271 4...
## 2 MULTIPOLYGON (((-87.60997 4...
## 3 MULTIPOLYGON (((-87.61633 4...
## 4 MULTIPOLYGON (((-87.63376 4...
## 5 MULTIPOLYGON (((-87.62064 4...
## 6 MULTIPOLYGON (((-87.63397 4...
2.7 Thematic Map Panel
To facilitate data discovery, we likely want to explore multiple maps at once. Here we’ll generate maps for multiple variables, and plot them as a map panel.
Can you think of more efficient ways to run this code? There are also other tmap
tricks to optimize this further, so enjoy your journey!
COVID <- tm_shape(Chi_Zipsf) + tm_fill("Case.Rate...Cumulative",
style="jenks", pal="BuPu", n=4, title = "COVID Rt") +
tm_layout(frame = F)
Senior <- tm_shape(Chi_Zipsf) + tm_fill("ovr65P",
style="jenks", pal="BuPu", n=4) +
tm_layout(frame = F)
NoHS <- tm_shape(Chi_Zipsf) + tm_fill("noHSP",
style="jenks", pal="BuPu", n=4) +
tm_layout(frame = F)
BlkP <- tm_shape(Chi_Zipsf) + tm_fill("blackP",
style="jenks", pal="BuPu", n=4) +
tm_layout(frame = F)
Latnx <- tm_shape(Chi_Zipsf) + tm_fill("hispP",
style="jenks", pal="BuPu", n=4) +
tm_layout(frame = F)
WhiP <- tm_shape(Chi_Zipsf) + tm_fill("whiteP",
style="jenks", pal="BuPu", n=4) +
tm_layout(frame = F)
tmap_arrange(COVID, Senior, NoHS, BlkP, Latnx, WhiP)
From the results, we see that cumulative COVID outcomes for one week in September 2020 seemed to have some geographic correlation with the Latinx/Hispanic community in Chicago. At the same time, low high school diploma rates are also concentrated in these areas, and there is some intersection with other variables considered. What are additional variables you could bring in to refine your approach? Perhaps percentage of essential workers; a different age group; internet access? What about linking in health outcomes like Asthma, Hypertension, and more at a similar scale?
In modern spatial epidemiology, associations must never be taken at face value. For example, we know that it is not “race” but “racism” that drives multiple health disparities – simply looking at a specific racial/ethnic group is not enough. Thus exploring multiple variables and nurturing a curiosity to understand these complex intersections will support knowledge discovery.
2.8 Write Data
We’re done! Well… not so fast. Let’s save the data so we don’t have to run the codebook again to access the data. Here, we’ll save as a geojson
file. This spatial format is more forgiving with long column names, which is a long-standing challenge with shapefiles. But sometimes it can be written oddly, so double-check the dimensions when reading in later.
We could also just write the data as a CSV
file which may be easier for future linking. For this, we use the st_drop_geometry()
function to remove the geometry data, so we’re just left with the attributes.
More Resources
For choropleth mapping in R:


Deadline looms for comments on plan to expand paid street parking in Cape Town | Cape Argus
Cape Town - Cape Town residents have been reminded to comment on a proposal to expand the areas subjected to paid street parking.
In a statement on Sunday, Felicity Purchase, mayco member for transport, told residents that the deadline to submit their comments and suggestions on the proposal to expand managed parking to more areas across Cape Town is Friday.
Sub-councils, business owners, and residents from CBDs and popular recreational areas like Camps Bay and Muizenberg have requested the City to expand managed parking to these areas due to the rapid increase in demand for on-street parking.
Purchase said: “Managed parking means the City charges a tariff for the use of an on-street parking bay; that the tariff is determined by the time the bay is being occupied and that those who park in CBDs are charged per 15 minutes.
“Visitors to recreational areas will be charged per three hours; that time restrictions are applicable; and the management times are from 8am to 5pm on weekdays and 8am to 1pm on Saturdays,” said Purchase.


blackjok3r
28-Nov-18 03:38 AM
kind of.... we tested a whole bunch of stuff at once with Ducks
there is a special mining mode where 10% fo the blcok reward mines to a multi sig address for making payments for notarisations
and we used it to select who could be a main net notary.
however... somwhere along the line MoMoM has broken again... it wasworking perfectly in the last stress test WITHOUT the CC mempool fix, and this is the only thing I can think of that is causing it ... so I might revisit that, and see if I can exempt notarisitions from it, like z-tx are exempt.
there is a way to esily check if a tx is a CC transaction., but in this case I dont think that will work, as a CC tx could depend on a normal tx, so I will need to use the code from the miner, that I made that gives notarisations top priotity, to detect notarisations and have them not be removed from the mempool. all starting to get very complicated tbh.... it should not ahve to work like this.


{% extends 'base.html' %} {% load compress %} {% block content %}
- move image around
- zoom with mouse scroll
Z - undo
X - redo
F - toggle polygon fill
A - enter adjust mode
D - enter draw mode
ESC - deselect any polygon and reset zoom
- add a new point
- close polygon
Q - snap to an existing line/boundary
W - snap to an existing point/corner
C - add a new point (same as left click)
- select existing polygon or point
- add new point to existing polygon
A - select next polygon to adjust
S - zoom in on selected polygon
DELETE or BACKSPACE - delete selected polygon (if no point under cursor)
DELETE or BACKSPACE - delete selected point (if point under cursor)


Ruby Classes: How and Why
Feb 29 2016
Arrays and Hashes are great for storing individual sets of data, but what if you need to store different versions of the same thing? For instance, if you want to track all of your favorite movies, you could create an individual hash for each containing the title, year, director, and lead actor, but that would take forever. Instead, you can use a Class to store this information, which gives you the access of a hash but also the flexibility to create your own.
To create a Class, you first need to define it with its initial arguments:
class Movie
initialize(title,year,director,actor)
end
end
This means that we have now created a new class called "Movie", and that whenever we log a new movie, it will have a title, year, director and actor. To make these data points generally available for use within the class, we will define them as instance variables:
class Movie
initialize(title,year,director,actor)
@title = title
@year = year
@director = director
@actor = actor
end
end
These 'instance variables', determined with the prefix @, are variables that can be at any time within the class, regardless of which movie (or instance of the class) you are using. These variables are not available outside of the class on their own, meaning other methods cannot find them, without a reader or writer method such as attr_reader, attr_writer, or attr_accessor, shown here:
class Movie
attr_reader :title, :year, :director, :actor
initialize(title,year,director,actor)
@title = title
@year = year
@director = director
@actor = actor
end
end
Now these instance variables can be called by other methods outside of the class.
You can also create methods specific to this class. For instance, if I wanted to create a shorthand way to view the title, year, and director for any given movie, I can create an instance method within the Movie class that will only work for my movies. But it will still work for any movie I put into my class. Check it out:
class Movie
attr_reader :title, :year, :director, :actor
initialize(title,year,director,actor)
@title = title
@year = year
@director = director
@actor = actor
end
def overview
puts "#{@title}, #{year}. Directed by {@director}. Stars<EMAIL_ADDRESS>end
end
This new overview method will now be available for any movies I put into this class, or any movies I enter in the future.
Classes are a great way to store multiple versions of the same thing, where you can define the precise methods and variables types you will need. You only have to set the class up once, and it will be very simple to add new instances over time without having to recode or define the same structure or methods over again.


Practices can increase referrals, reduce the number of phone calls and see more patients. Implement the pingmd business model that suits your practice and stand out.
Deliver Timely Care
ACOs and Hospitals need to coordinate care with their care teams, increasing their capacity to engage with more patients and deliver timely care decisions to reduce high cost negative events. Pingmd optimizes your network and teams.
Engage with patients via e-visit
High quality care is made convenient on pingmd. Involve your care network, your family and see how easy it is to keep track of your visit instructions and reach out to your care team in a way that they are more able to serve you when you really need it.
Expand Your Care Network
High quality care is made convenient on pingmd. Involve your care network, your family and see how easy it is to keep track of your visit instructions and reach out to your care team in a way that they are more able to serve you when you really need it.
Contact Us
Any questions regarding plans, our support team can help lorem ipsum dolor sit amet quidam mutat qui te.