## Code to generate and readily evaluated nodes and weights

To accompany the paper "Likelihood Approximation by Numerical Integration on Sparse Grids"
by Florian Heiss & Viktor Winschel

This website offers three options to use quadrature on sparse grids:

Nov. 12, 2007:
The first versions are online. If you find bugs or have comments or suggestions, let us know!

## General remarks

To approximate a D-variate integral of the form
Ωg(x)f(x)dx
where Ω⊆RD and f(x) is a weight function, a quadrature rule prescribes a RxD matrix of nodes X and a Rx1 vector of weights w depending on Ω and f. The approximation is calculated as
g(X)'w.

In all three versions provided here, four quadrature rules for D-dimensional numerical integration are implemented:

• Rules for the unweighted integration: Ω=[0,1]D, f(x)=1
• GQU: Univariate Gaussian quadrature rules as basis
• KPU: Univariate nested quadrature rules as basis - delayed Kronrod-Patterson rules, see Knut Petras (2003): "Smolyak cubature of given polynomial degree with few nodes for increasing dimension." Numerische Mathematik 93, 729-753.
• Rules for the integration over Ω=RD with Gaussian weight f(x)=exp(-x'x/2)/sqrt(2π)
• GQN: Univariate Gaussian quadrature rules as basis
• KPN: Univariate nested quadrature rules as basis, see A. Genz and B. D. Keister (1996): "Fully symmetric interpolatory rules for multiple integrals over infinite regions with Gaussian weight." Journal of Computational and Applied Mathematics 71, 299-309.

As a general rule, the versions with the nested univariate rules are more efficient in the sense that they require fewer nodes to achieve the same level of polynomial exactness. The Matlab code and the Stata code allow the user to use own one-dimensional quadrature rules.

The accuracy level in each of the implementations is denoted as k in the paper. Each rule with accuracy level k integrates complete polynomials of total order 2k-1 exactly.

## Matlab code

The Matlab code can be used to extend any univariate quadrature rule to multiple dimensions. The four rules discussed above are included in the code.

To install the package, perform the following steps:

• Make sure Matlab will find the file: Put it into the current working directory, or somewhere into the Matlab search path. You can view and modify it in the menu, see `File``Set Path...`

The syntax is the following:
`[n,w]=nwspgr(type,dim,acc,sym)`

• `type` is a string choosing the underlying univariate quadrature rule. It is either `'GQU'`, `'KPU'`, `'GQN'`, or `'KPN'` for the built-in rules, or the name of a function for user-supplied univariate rules. In the latter case, the user must make sure that this function is known to Matlab and that it
• accepts one scalar argument as the accuracy level
• returns a column vector of nodes and a column vector of weights.
for example this function could be defined as `function [n w] = myquadrule(acc)`. In this case, the argument `type` would be the string `'myquadrule'`.
• `dim` is the number of dimensions.
• `acc` is the accuracy level.
• `sym` is optional and only interpreted if a user-supplied univariate quadrature rule is used. In this case, if `sym` is used and a value other than 0 is supplied, this has two effects: (1) the code runs much faster with high `dim` and/or `acc` and (2) the result will only be valid if the underlying univariate rule is "symmetric". This is the case if
• it uses one node for accuracy level 1, call it m
• in higher accuracy levels, the nodes cosist of m and/or a set of other nodes. For each of these other nodes m+d, there is also a node m-d with the same weight.
• The results will be returned as a (Rx`dim`) matrix of nodes `n` and a (Rx1) vector of corresponding weights `w`, where R depends on `dim` and `acc`, see the examples below.

An example of the program in action is in nwSpGr_demo.m.

## Stata code

The Stata code will run in Stata versions 9 and above. There are three ways the nodes and weights can be obtained:

• As a Mata matrix
• As a Stata matrix
• As variables

To install the package, perform the following steps:

• Download stata.zip and unpack its contents to a place where Stata will find the files. A natural place is the personal ado directory. If you don't know where it is, Stata's command `sysdir` will tell you.
• In order to re-index the installed functions, either close and reopen Stata or use the command `mata: mata mlib index`.

The syntax is the following:

• For generating nodes and weights as Stata matrices and/or variables, use the command
```nwspgr, type(string) dimensions(integer) accuracy(integer) [ matnodes(string) matweights(string) varnodes(string) varweights(string) ] ```
Example: `nwspgr, type(KPN) dim(2) acc(4) varn(x) varw(w)`
The first three options determine what is calculated:
• `type(string)` chooses the underlying univariate quadrature rule, where `string` is either `GQU`, `KPU`, `GQN`, or `KPN`.
• `dim(integer)` determines the number of dimensions.
• `acc(integer)` determines the accuracy level.
The other options determine how the results are returned:
• `matn(string)` requests that the matrix of nodes is returned as a matrix named `string`.
• `matw(string)` requests that the vector of weights is returned as a matrix named `string`.
• `varn(string)` requests that the matrix of nodes is returned as variables named `string1`,...,`stringD` .
• `varw(string)` requests that the vector of weights is returned as a variable named `string`.
Not specifying any of those options is equivalent to specifying `matn(nodes) matw(weights)`.
• For generating a matrix of nodes and weights from within Stata's matrix language Mata, use the Mata function
`nwspgr(type,dim,acc,sym)`
Example: `mata: nw = nwspgr("KPN",2,4)`
• `type` chooses the underlying univariate quadrature rule. It can be one of the strings `"GQU"`, `"KPU"`, `"GQN"`, or `"KPN"` if the built-in rules are used. If the user wants to use an own underlying univariate quadrature rule, `type` has to be a pointer to a function. This function
• accepts one scalar argument as the accuracy level
• returns a matrix with two columns, where the first column represents nodes and the second column represents weights.
for example this function could be declared as `real matrix myquadrule(real scalar acc)`. In this case, the argument `type` would be the pointer `&myquadrule()`.
• `dim` is an integer determining the number of dimensions.
• `acc` is an integer determining the accuracy level.
• `sym` is optional and only interpreted if a user-supplied univariate quadrature rule is used. In this case, if `sym` is used and a value other than 0 is supplied, this has two effects: (1) the code runs much faster with high `dim` and/or `acc` and (2) the result will only be valid if the underlying univariate rule is "symmetric". This is the case if
• it uses one node for accuracy level 1, call it m
• in higher accuracy levels, the nodes cosist of m and/or a set of other nodes. For each of these other nodes m+d, there is also a node m-d with the same weight.
The result is a Rx`dim+1` matrix where the last column represents the weights.

One example for each of the types of use is here:

The source code of the compiled Mata library `lnwspgr.mlib` contained in the installation ZIP file is not necessary to use the program. For the interest users, here it is: build_nwspgr.do

## Readily evaluated nodes and weights

For a variety of combinations of the number of dimensions and the accuracy level, readily generated nodes and weights can be downloaded as comma-seperated ASCII files. The file names are
TYPE_dDIM_lLEV.asc
where TYPE denotes the quadrature rule and is GQU, KPU, GQN, or KPN, DIM is the number of dimensions with 1≤DIM≤20, and LEV is the accuracy level with 1≤LEV≤25. The ASCII files for DIM dimensions have DIM+1 columns, with the last column representing the integration weight.

The following tables show the number of nodes required by the combinations of TYPE, DIM, and LEV. Note that since the number of nodes increases with the number of dimensions, rules for higher dimensions are only implemented for lower accuracy levels. The files are provides as ZIPped archives - one for each TYPE and for each combination of TYPE and DIM.

### Rule = KPU

ZIP file with all files: KPU.zip
Table of number of nodes with links to separate files for each number of dimensions:
Accuracy
Dim. 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25
1 1 3 3 7 7 7 15 15 15 15 15 15 31 31 31 31 31 31 31 31 31 31 31 31 63
2 1 5 9 17 33 33 65 97 97 161 161 161 257 321 321 449 449 449 705 705 705 705 705 705 1025
3 1 7 19 39 87 135 207 399 495 751 1135 1135 1759 2335 2527 3679 4447 4447 6495 8031 8031 11103 11103 11103 15039
4 1 9 33 81 193 385 641 1217 1985 2881 4929 6465 8705 13697 16001 22401 31617 34689 47489 63873
5 1 11 51 151 391 903 1743 3343 6223 10063 17103 27343 38303 60703
6 1 13 73 257 737 1889 4161 8481 16929 30689 53729 93665
7 1 15 99 407 1303 3655 8975 19855 42031 83247
8 1 17 129 609 2177 6657 17921 43137 97153
9 1 19 163 871 3463 11527 33679 87823
10 1 21 201 1201 5281 19105 60225
11 1 23 243 1607 7767 30471
12 1 25 289 2097 11073 46977
13 1 27 339 2679 15367 70279
14 1 29 393 3361 20833 102369
15 1 31 451 4151 27671
16 1 33 513 5057 36097
17 1 35 579 6087 46343
18 1 37 649 7249 58657
19 1 39 723 8551 73303
20 1 41 801 10001 90561

### Rule = KPN

ZIP file with all files: KPN.zip
Table of number of nodes with links to separate files for each number of dimensions:
Accuracy
Dim. 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25
1 1 3 3 7 9 9 9 9 17 19 19 19 19 19 19 31 33 35 35 35 35 35 35 35 35
2 1 5 9 17 37 45 61 77 97 133 141 205 253 261 261 285 401 445 553 617 641 649 649 841 921
3 1 7 19 39 93 165 237 381 513 703 919 1183 1719 2031 2463 2979 3513 4191 4731 6315 7539 8523 9267 10179 12203
4 1 9 33 81 201 441 761 1305 2129 3065 4489 6185 8745 12057 15321 20681 25985 32025 39233 48321
5 1 11 51 151 401 993 2033 3793 6913 11323 17643 27003 39403 57563
6 1 13 73 257 749 2021 4725 9765 19281 35357 59957 98837
7 1 15 99 407 1317 3837 9941 22725 48401 96967
8 1 17 129 609 2193 6897 19441 48689 111841
9 1 19 163 871 3481 11833 35929 97561
10 1 21 201 1201 5301 19485 63405
11 1 23 243 1607 7789 30933
12 1 25 289 2097 11097 47529
13 1 28 339 2679 15393 70929
14 1 29 393 3361 20861 103125
15 1 31 451 4151 27701
16 1 33 513 5057 36129
17 1 35 579 6087 46377
18 1 37 649 7249 58693
19 1 39 723 8551 73341
20 1 41 801 10001 90601

### Rule = GQU

ZIP file with all files: GQU.zip
Table of number of nodes with links to separate files for each number of dimensions:
Accuracy
Dim. 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25
1 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25
2 1 5 13 29 53 89 137 201 281 381 501 645 813 1009 1233 1489 1777 2101 2461 2861 3301 3785 4313 4889 5513
3 1 7 25 69 165 351 681 1233 2097 3407 5297 7973 11621 16535 22961 31297 41857 55159 71593
4 1 9 41 138 385 953 2145 4481 8785 16345 29033 49577 81713
5 1 11 61 241 781 2203 5593 13073 28553 58923 115813
6 1 13 85 389 1433 4541 12841 33193 79729
7 1 15 113 589 2437 8583 26769 75841
8 1 17 145 849 3905 15153 51713
9 1 19 181 1177 5965 25315 93913
10 1 21 221 1581 8761 40405 162025
11 1 23 265 2069 12453 62063
12 1 25 313 2649 17217 92265
13 1 27 365 3329 23245 133355
14 1 29 421 4117 30745
15 1 31 481 5021 39941
16 1 33 545 6049 51073
17 1 35 613 7209 64397
18 1 37 685 8509 80185
19 1 39 761 9957 98725
20 1 41 841 11561 120321

### Rule = GQN

ZIP file with all files: GQN.zip
Table of number of nodes with links to separate files for each number of dimensions:
Accuracy
Dim. 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25
1 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25
2 1 5 13 29 53 89 137 201 281 381 501 645 813 1009 1233 1489 1777 2101 2461 2861 3301 3785 4313 4889 5513
3 1 7 25 69 165 351 681 1233 2097 3407 5297 7973 11621 16535 22961 31297 41857 55159 71593
4 1 9 41 137 385 953 2145 4481 8785 16345 29033 49577 81713
5 1 11 61 241 781 2203 5593 13073 28553 58923 115813
6 1 13 85 389 1433 4541 12841 33193 79729
7 1 15 113 589 2437 8583 26769 75841
8 1 17 145 849 3905 15153 51713
9 1 19 181 1177 5965 25315 93913
10 1 21 221 1581 8761 40405 162025
11 1 23 265 2069 12453 62063
12 1 25 313 2649 17217 92265
13 1 27 365 3329 23245 133355
14 1 29 421 4117 30745
15 1 31 481 5021 39941
16 1 33 545 6049 51073
17 1 35 613 7209 64397
18 1 37 685 8509 80185
19 1 39 761 9957 98725
20 1 41 841 11561 120321

### Contact:

 Florian Heiss Viktor Winschel University of Munich University of Mannheim Department of Economics Department of Economics f.heiss@lmu.de winschel@rumms.uni-mannheim.de

Updated Nov 12, 2007