text:
#==============================================================================#
# CloudSearchDomain.jl
#
# This file is generated from:
# https://github.com/aws/aws-sdk-js/blob/master/apis/cloudsearchdomain-2013-01-01.normal.json
#==============================================================================#

__precompile__()

module CloudSearchDomain

using AWSCore

"""
    using AWSSDK.CloudSearchDomain.search
    search([::AWSConfig], arguments::Dict)
    search([::AWSConfig]; q=, <keyword arguments>)

    using AWSCore.Services.cloudsearchdomain
    cloudsearchdomain([::AWSConfig], "GET", "/2013-01-01/search?format=sdk&pretty=true", arguments::Dict)
    cloudsearchdomain([::AWSConfig], "GET", "/2013-01-01/search?format=sdk&pretty=true", q=, <keyword arguments>)

# Search Operation

Retrieves a list of documents that match the specified search criteria. How you specify the search criteria depends on which query parser you use. Amazon CloudSearch supports four query parsers:

*   `simple`: search all `text` and `text-array` fields for the specified string. Search for phrases, individual terms, and prefixes.
*   `structured`: search specific fields, construct compound queries using Boolean operators, and use advanced features such as term boosting and proximity searching.
*   `lucene`: specify search criteria using the Apache Lucene query parser syntax.
*   `dismax`: specify search criteria using the simplified subset of the Apache Lucene query parser syntax defined by the DisMax query parser.

For more information, see [Searching Your Data](http://docs.aws.amazon.com/cloudsearch/latest/developerguide/searching.html) in the *Amazon CloudSearch Developer Guide*.

The endpoint for submitting `Search` requests is domain-specific. You submit search requests to a domain's search endpoint. To get the search endpoint for your domain, use the Amazon CloudSearch configuration service `DescribeDomains` action. A domain's endpoints are also displayed on the domain dashboard in the Amazon CloudSearch console.

# Arguments

## `cursor = ::String`
Retrieves a cursor value you can use to page through large result sets. Use the `size` parameter to control the number of hits to include in each response. You can specify either the `cursor` or `start` parameter in a request; they are mutually exclusive. To get the first cursor, set the cursor value to `initial`. In subsequent requests, specify the cursor value returned in the hits section of the response.

For more information, see [Paginating Results](http://docs.aws.amazon.com/cloudsearch/latest/developerguide/paginating-results.html) in the *Amazon CloudSearch Developer Guide*.

## `expr = ::String`
Defines one or more numeric expressions that can be used to sort results or specify search or filter criteria. You can also specify expressions as return fields.

You specify the expressions in JSON using the form `{"EXPRESSIONNAME":"EXPRESSION"}`. You can define and use multiple expressions in a search request. For example:

`{"expression1":"_score*rating", "expression2":"(1/rank)*year"}`

For information about the variables, operators, and functions you can use in expressions, see [Writing Expressions](http://docs.aws.amazon.com/cloudsearch/latest/developerguide/configuring-expressions.html#writing-expressions) in the *Amazon CloudSearch Developer Guide*.

## `facet = ::String`
Specifies one or more fields for which to get facet information, and options that control how the facet information is returned. Each specified field must be facet-enabled in the domain configuration. The fields and options are specified in JSON using the form `{"FIELD":{"OPTION":VALUE,"OPTION:"STRING"},"FIELD":{"OPTION":VALUE,"OPTION":"STRING"}}`.

You can specify the following faceting options:

*   `buckets` specifies an array of the facet values or ranges to count. Ranges are specified using the same syntax that you use to search for a range of values. For more information, see [Searching for a Range of Values](http://docs.aws.amazon.com/cloudsearch/latest/developerguide/searching-ranges.html) in the *Amazon CloudSearch Developer Guide*. Buckets are returned in the order they are specified in the request. The `sort` and `size` options are not valid if you specify `buckets`.

*   `size` specifies the maximum number of facets to include in the results. By default, Amazon CloudSearch returns counts for the top 10. The `size` parameter is only valid when you specify the `sort` option; it cannot be used in conjunction with `buckets`.

*   `sort` specifies how you want to sort the facets in the results: `bucket` or `count`. Specify `bucket` to sort alphabetically or numerically by facet value (in ascending order). Specify `count` to sort by the facet counts computed for each facet value (in descending order). To retrieve facet counts for particular values or ranges of values, use the `buckets` option instead of `sort`.

If no facet options are specified, facet counts are computed for all field values, the facets are sorted by facet count, and the top 10 facets are returned in the results.

To count particular buckets of values, use the `buckets` option. For example, the following request uses the `buckets` option to calculate and return facet counts by decade.

`{"year":{"buckets":["[1970,1979]","[1980,1989]","[1990,1999]","[2000,2009]","[2010,}"]}}`

To sort facets by facet count, use the `count` option. For example, the following request sets the `sort` option to `count` to sort the facet values by facet count, with the facet values that have the most matching documents listed first. Setting the `size` option to 3 returns only the top three facet values.

`{"year":{"sort":"count","size":3}}`

To sort the facets by value, use the `bucket` option. For example, the following request sets the `sort` option to `bucket` to sort the facet values numerically by year, with earliest year listed first.

`{"year":{"sort":"bucket"}}`

For more information, see [Getting and Using Facet Information](http://docs.aws.amazon.com/cloudsearch/latest/developerguide/faceting.html) in the *Amazon CloudSearch Developer Guide*.

## `fq = ::String`
Specifies a structured query that filters the results of a search without affecting how the results are scored and sorted. You use `filterQuery` in conjunction with the `query` parameter to filter the documents that match the constraints specified in the `query` parameter. Specifying a filter controls only which matching documents are included in the results, it has no effect on how they are scored and sorted. The `filterQuery` parameter supports the full structured query syntax.

For more information about using filters, see [Filtering Matching Documents](http://docs.aws.amazon.com/cloudsearch/latest/developerguide/filtering-results.html) in the *Amazon CloudSearch Developer Guide*.

## `highlight = ::String`
Retrieves highlights for matches in the specified `text` or `text-array` fields. Each specified field must be highlight enabled in the domain configuration. The fields and options are specified in JSON using the form `{"FIELD":{"OPTION":VALUE,"OPTION:"STRING"},"FIELD":{"OPTION":VALUE,"OPTION":"STRING"}}`.

You can specify the following highlight options:

*   `format`: specifies the format of the data in the text field: `text` or `html`. When data is returned as HTML, all non-alphanumeric characters are encoded. The default is `html`.
*   `max_phrases`: specifies the maximum number of occurrences of the search term(s) you want to highlight. By default, the first occurrence is highlighted.
*   `pre_tag`: specifies the string to prepend to an occurrence of a search term. The default for HTML highlights is `<em>`. The default for text highlights is `*`.
*   `post_tag`: specifies the string to append to an occurrence of a search term. The default for HTML highlights is `</em>`. The default for text highlights is `*`.

If no highlight options are specified for a field, the returned field text is treated as HTML and the first match is highlighted with emphasis tags: `<em>search-term</em>`.

For example, the following request retrieves highlights for the `actors` and `title` fields.

`{ "actors": {}, "title": {"format": "text","max_phrases": 2,"pre_tag": "**","post_tag": "**"} }`

## `partial = ::Bool`
Enables partial results to be returned if one or more index partitions are unavailable. When your search index is partitioned across multiple search instances, by default Amazon CloudSearch only returns results if every partition can be queried. This means that the failure of a single search instance can result in 5xx (internal server) errors. When you enable partial results, Amazon CloudSearch returns whatever results are available and includes the percentage of documents searched in the search results (percent-searched). This enables you to more gracefully degrade your users' search experience. For example, rather than displaying no results, you could display the partial results and a message indicating that the results might be incomplete due to a temporary system outage.

## `q = ::String` -- *Required*
Specifies the search criteria for the request. How you specify the search criteria depends on the query parser used for the request and the parser options specified in the `queryOptions` parameter. By default, the `simple` query parser is used to process requests. To use the `structured`, `lucene`, or `dismax` query parser, you must also specify the `queryParser` parameter.

For more information about specifying search criteria, see [Searching Your Data](http://docs.aws.amazon.com/cloudsearch/latest/developerguide/searching.html) in the *Amazon CloudSearch Developer Guide*.

## `q.options = ::String`
Configures options for the query parser specified in the `queryParser` parameter. You specify the options in JSON using the following form `{"OPTION1":"VALUE1","OPTION2":VALUE2"..."OPTIONN":"VALUEN"}.`

The options you can configure vary according to which parser you use:

*   `defaultOperator`: The default operator used to combine individual terms in the search string. For example: `defaultOperator: 'or'`. For the `dismax` parser, you specify a percentage that represents the percentage of terms in the search string (rounded down) that must match, rather than a default operator. A value of `0%` is the equivalent to OR, and a value of `100%` is equivalent to AND. The percentage must be specified as a value in the range 0-100 followed by the percent (%) symbol. For example, `defaultOperator: 50%`. Valid values: `and`, `or`, a percentage in the range 0%-100% (`dismax`). Default: `and` (`simple`, `structured`, `lucene`) or `100` (`dismax`). Valid for: `simple`, `structured`, `lucene`, and `dismax`.
*   `fields`: An array of the fields to search when no fields are specified in a search. If no fields are specified in a search and this option is not specified, all text and text-array fields are searched. You can specify a weight for each field to control the relative importance of each field when Amazon CloudSearch calculates relevance scores. To specify a field weight, append a caret (`^`) symbol and the weight to the field name. For example, to boost the importance of the `title` field over the `description` field you could specify: `"fields":["title^5","description"]`. Valid values: The name of any configured field and an optional numeric value greater than zero. Default: All `text` and `text-array` fields. Valid for: `simple`, `structured`, `lucene`, and `dismax`.
*   `operators`: An array of the operators or special characters you want to disable for the simple query parser. If you disable the `and`, `or`, or `not` operators, the corresponding operators (`+`, `|`, `-`) have no special meaning and are dropped from the search string. Similarly, disabling `prefix` disables the wildcard operator (`*`) and disabling `phrase` disables the ability to search for phrases by enclosing phrases in double quotes. Disabling precedence disables the ability to control order of precedence using parentheses. Disabling `near` disables the ability to use the ~ operator to perform a sloppy phrase search. Disabling the `fuzzy` operator disables the ability to use the ~ operator to perform a fuzzy search. `escape` disables the ability to use a backslash (`\\`) to escape special characters within the search string. Disabling whitespace is an advanced option that prevents the parser from tokenizing on whitespace, which can be useful for Vietnamese. (It prevents Vietnamese words from being split incorrectly.) For example, you could disable all operators other than the phrase operator to support just simple term and phrase queries: `"operators":["and","not","or", "prefix"]`. Valid values: `and`, `escape`, `fuzzy`, `near`, `not`, `or`, `phrase`, `precedence`, `prefix`, `whitespace`. Default: All operators and special characters are enabled. Valid for: `simple`.
*   `phraseFields`: An array of the `text` or `text-array` fields you want to use for phrase searches. When the terms in the search string appear in close proximity within a field, the field scores higher. You can specify a weight for each field to boost that score. The `phraseSlop` option controls how much the matches can deviate from the search string and still be boosted. To specify a field weight, append a caret (`^`) symbol and the weight to the field name. For example, to boost phrase matches in the `title` field over the `abstract` field, you could specify: `"phraseFields":["title^3", "plot"]` Valid values: The name of any `text` or `text-array` field and an optional numeric value greater than zero. Default: No fields. If you don't specify any fields with `phraseFields`, proximity scoring is disabled even if `phraseSlop` is specified. Valid for: `dismax`.
*   `phraseSlop`: An integer value that specifies how much matches can deviate from the search phrase and still be boosted according to the weights specified in the `phraseFields` option; for example, `phraseSlop: 2`. You must also specify `phraseFields` to enable proximity scoring. Valid values: positive integers. Default: 0. Valid for: `dismax`.
*   `explicitPhraseSlop`: An integer value that specifies how much a match can deviate from the search phrase when the phrase is enclosed in double quotes in the search string. (Phrases that exceed this proximity distance are not considered a match.) For example, to specify a slop of three for dismax phrase queries, you would specify `"explicitPhraseSlop":3`. Valid values: positive integers. Default: 0. Valid for: `dismax`.
*   `tieBreaker`: When a term in the search string is found in a document's field, a score is calculated for that field based on how common the word is in that field compared to other documents. If the term occurs in multiple fields within a document, by default only the highest scoring field contributes to the document's overall score. You can specify a `tieBreaker` value to enable the matches in lower-scoring fields to contribute to the document's score. That way, if two documents have the same max field score for a particular term, the score for the document that has matches in more fields will be higher. The formula for calculating the score with a tieBreaker is `(max field score) + (tieBreaker) * (sum of the scores for the rest of the matching fields)`. Set `tieBreaker` to 0 to disregard all but the highest scoring field (pure max): `"tieBreaker":0`. Set to 1 to sum the scores from all fields (pure sum): `"tieBreaker":1`. Valid values: 0.0 to 1.0. Default: 0.0. Valid for: `dismax`.

## `q.parser = "simple", "structured", "lucene" or "dismax"`
Specifies which query parser to use to process the request. If `queryParser` is not specified, Amazon CloudSearch uses the `simple` query parser.

Amazon CloudSearch supports four query parsers:

*   `simple`: perform simple searches of `text` and `text-array` fields. By default, the `simple` query parser searches all `text` and `text-array` fields. You can specify which fields to search by with the `queryOptions` parameter. If you prefix a search term with a plus sign (+) documents must contain the term to be considered a match. (This is the default, unless you configure the default operator with the `queryOptions` parameter.) You can use the `-` (NOT), `|` (OR), and `*` (wildcard) operators to exclude particular terms, find results that match any of the specified terms, or search for a prefix. To search for a phrase rather than individual terms, enclose the phrase in double quotes. For more information, see [Searching for Text](http://docs.aws.amazon.com/cloudsearch/latest/developerguide/searching-text.html) in the *Amazon CloudSearch Developer Guide*.
*   `structured`: perform advanced searches by combining multiple expressions to define the search criteria. You can also search within particular fields, search for values and ranges of values, and use advanced options such as term boosting, `matchall`, and `near`. For more information, see [Constructing Compound Queries](http://docs.aws.amazon.com/cloudsearch/latest/developerguide/searching-compound-queries.html) in the *Amazon CloudSearch Developer Guide*.
*   `lucene`: search using the Apache Lucene query parser syntax. For more information, see [Apache Lucene Query Parser Syntax](http://lucene.apache.org/core/4_6_0/queryparser/org/apache/lucene/queryparser/classic/package-summary.html#package_description).
*   `dismax`: search using the simplified subset of the Apache Lucene query parser syntax defined by the DisMax query parser. For more information, see [DisMax Query Parser Syntax](http://wiki.apache.org/solr/DisMaxQParserPlugin#Query_Syntax).

## `return = ::String`
Specifies the field and expression values to include in the response. Multiple fields or expressions are specified as a comma-separated list. By default, a search response includes all return enabled fields (`_all_fields`). To return only the document IDs for the matching documents, specify `_no_fields`. To retrieve the relevance score calculated for each document, specify `_score`.

## `size = ::Int`
Specifies the maximum number of search hits to include in the response.

## `sort = ::String`
Specifies the fields or custom expressions to use to sort the search results. Multiple fields or expressions are specified as a comma-separated list. You must specify the sort direction (`asc` or `desc`) for each field; for example, `year desc,title asc`. To use a field to sort results, the field must be sort-enabled in the domain configuration. Array type fields cannot be used for sorting. If no `sort` parameter is specified, results are sorted by their default relevance scores in descending order: `_score desc`. You can also sort by document ID (`_id asc`) and version (`_version desc`).

For more information, see [Sorting Results](http://docs.aws.amazon.com/cloudsearch/latest/developerguide/sorting-results.html) in the *Amazon CloudSearch Developer Guide*.

## `start = ::Int`
Specifies the offset of the first search hit you want to return. Note that the result set is zero-based; the first result is at index 0. You can specify either the `start` or `cursor` parameter in a request, they are mutually exclusive.

For more information, see [Paginating Results](http://docs.aws.amazon.com/cloudsearch/latest/developerguide/paginating-results.html) in the *Amazon CloudSearch Developer Guide*.

## `stats = ::String`
Specifies one or more fields for which to get statistics information. Each specified field must be facet-enabled in the domain configuration. The fields are specified in JSON using the form:

`{"FIELD-A":{},"FIELD-B":{}}`

There are currently no options supported for statistics.

# Returns

`SearchResponse`

# Exceptions

`SearchException`.

See also: [AWS API Documentation](https://docs.aws.amazon.com/goto/WebAPI/cloudsearchdomain-2013-01-01/Search)
"""
@inline search(aws::AWSConfig=default_aws_config(); args...) = search(aws, args)

@inline search(aws::AWSConfig, args) = AWSCore.Services.cloudsearchdomain(aws, "GET", "/2013-01-01/search?format=sdk&pretty=true", args)

@inline search(args) = search(default_aws_config(), args)

"""
    using AWSSDK.CloudSearchDomain.suggest
    suggest([::AWSConfig], arguments::Dict)
    suggest([::AWSConfig]; q=, suggester=, <keyword arguments>)

    using AWSCore.Services.cloudsearchdomain
    cloudsearchdomain([::AWSConfig], "GET", "/2013-01-01/suggest?format=sdk&pretty=true", arguments::Dict)
    cloudsearchdomain([::AWSConfig], "GET", "/2013-01-01/suggest?format=sdk&pretty=true", q=, suggester=, <keyword arguments>)

# Suggest Operation

Retrieves autocomplete suggestions for a partial query string. You can use suggestions enable you to display likely matches before users finish typing. In Amazon CloudSearch, suggestions are based on the contents of a particular text field. When you request suggestions, Amazon CloudSearch finds all of the documents whose values in the suggester field start with the specified query string. The beginning of the field must match the query string to be considered a match.

For more information about configuring suggesters and retrieving suggestions, see [Getting Suggestions](http://docs.aws.amazon.com/cloudsearch/latest/developerguide/getting-suggestions.html) in the *Amazon CloudSearch Developer Guide*.

The endpoint for submitting `Suggest` requests is domain-specific. You submit suggest requests to a domain's search endpoint. To get the search endpoint for your domain, use the Amazon CloudSearch configuration service `DescribeDomains` action. A domain's endpoints are also displayed on the domain dashboard in the Amazon CloudSearch console.

# Arguments

## `q = ::String` -- *Required*
Specifies the string for which you want to get suggestions.

## `suggester = ::String` -- *Required*
Specifies the name of the suggester to use to find suggested matches.

## `size = ::Int`
Specifies the maximum number of suggestions to return.

# Returns

`SuggestResponse`

# Exceptions

`SearchException`.

See also: [AWS API Documentation](https://docs.aws.amazon.com/goto/WebAPI/cloudsearchdomain-2013-01-01/Suggest)
"""
@inline suggest(aws::AWSConfig=default_aws_config(); args...) = suggest(aws, args)

@inline suggest(aws::AWSConfig, args) = AWSCore.Services.cloudsearchdomain(aws, "GET", "/2013-01-01/suggest?format=sdk&pretty=true", args)

@inline suggest(args) = suggest(default_aws_config(), args)

"""
    using AWSSDK.CloudSearchDomain.upload_documents
    upload_documents([::AWSConfig], arguments::Dict)
    upload_documents([::AWSConfig]; documents=, *header:* Content-Type=)

    using AWSCore.Services.cloudsearchdomain
    cloudsearchdomain([::AWSConfig], "POST", "/2013-01-01/documents/batch?format=sdk", arguments::Dict)
    cloudsearchdomain([::AWSConfig], "POST", "/2013-01-01/documents/batch?format=sdk", documents=, *header:* Content-Type=)

# UploadDocuments Operation

Posts a batch of documents to a search domain for indexing. A document batch is a collection of add and delete operations that represent the documents you want to add, update, or delete from your domain. Batches can be described in either JSON or XML. Each item that you want Amazon CloudSearch to return as a search result (such as a product) is represented as a document. Every document has a unique ID and one or more fields that contain the data that you want to search and return in results. Individual documents cannot contain more than 1 MB of data. The entire batch cannot exceed 5 MB. To get the best possible upload performance, group add and delete operations in batches that are close the 5 MB limit. Submitting a large volume of single-document batches can overload a domain's document service.

The endpoint for submitting `UploadDocuments` requests is domain-specific. To get the document endpoint for your domain, use the Amazon CloudSearch configuration service `DescribeDomains` action. A domain's endpoints are also displayed on the domain dashboard in the Amazon CloudSearch console.

For more information about formatting your data for Amazon CloudSearch, see [Preparing Your Data](http://docs.aws.amazon.com/cloudsearch/latest/developerguide/preparing-data.html) in the *Amazon CloudSearch Developer Guide*. For more information about uploading data for indexing, see [Uploading Data](http://docs.aws.amazon.com/cloudsearch/latest/developerguide/uploading-data.html) in the *Amazon CloudSearch Developer Guide*.

# Arguments

## `documents = blob` -- *Required*
A batch of documents formatted in JSON or HTML.

## `*header:* Content-Type = "application/json" or "application/xml"` -- *Required*
The format of the batch you are uploading. Amazon CloudSearch supports two document batch formats:

*   application/json
*   application/xml

# Returns

`UploadDocumentsResponse`

# Exceptions

`DocumentServiceException`.

See also: [AWS API Documentation](https://docs.aws.amazon.com/goto/WebAPI/cloudsearchdomain-2013-01-01/UploadDocuments)
"""
@inline upload_documents(aws::AWSConfig=default_aws_config(); args...) = upload_documents(aws, args)

@inline upload_documents(aws::AWSConfig, args) = AWSCore.Services.cloudsearchdomain(aws, "POST", "/2013-01-01/documents/batch?format=sdk", args)

@inline upload_documents(args) = upload_documents(default_aws_config(), args)

end # module CloudSearchDomain

#==============================================================================#
# End of file
#==============================================================================#


text:
push!(LOAD_PATH,"../src/")
using Documenter, ImageBinarization
makedocs(sitename="Documentation",
    format = Documenter.HTML(
        prettyurls = get(ENV, "CI", nothing) == "true"
    ),
    pages = [
        "Home" => "index.md",
        "Function Reference" => "reference.md"
    ]
)
deploydocs(repo = "github.com/zygmuntszpak/ImageBinarization.jl.git")

ext: jl


text:
using Test, Econometrics, CSV, RDatasets, LinearAlgebra

# Linear Models
@testset "Linear Models" begin
	@testset "Balanced Panel Data" begin
    # Balanced Panel Data
    data = dataset("Ecdat", "Crime")
    # Pooling
    model = fit(EconometricModel,
                @formula(CRMRTE ~ PrbConv + AvgSen + PrbPris),
                data)
    @test nullloglikelihood(model) ≈ 1633.305 rtol = 1e-3
    @test loglikelihood(model) ≈ 1643.304 rtol = 1e-3
    @test dof(model) == 5
    # This corresponds to BIC from R's lm, Stata uses rank(vcov(model))
    @test bic(model) ≈ -3254.379 rtol = 1e-3
    @test coef(model) ≈ [0.0186413, -0.0011647, 0.0002362, 0.0273394] rtol = 1e-5
    @test vcov(model) ≈ [  1.858177e-05 -1.620480e-07 -6.464872e-07 -2.860941e-05
                          -1.620480e-07  1.781365e-07 -1.715589e-09  1.286486e-07
                          -6.464872e-07 -1.715589e-09  7.193971e-08  8.182148e-09
                          -2.860941e-05  1.286486e-07  8.182148e-09  6.685387e-05
                        ] rtol = 1e-6
    @test confint(model) ≈ [ 0.0101761950  0.0271063941
                            -0.0019935581 -0.0003358998
                            -0.0002905263  0.0007628970
                             0.0112828681  0.0433959412
                           ] rtol = 1e-6
    # Between
    model = fit(EconometricModel,
                @formula(CRMRTE ~ PrbConv + AvgSen + PrbPris + between(County)),
                data)
    @test nullloglikelihood(model) ≈ 239.5638 rtol = 1e-3
    @test loglikelihood(model) ≈ 244.4479 rtol = 1e-3
    @test dof(model) == 5
    @test coef(model) ≈ [-0.00464339, -0.00354113, 0.00084805, 0.07302994] rtol = 1e-6
    @test vcov(model) ≈ [ 3.480883e-04 -5.985740e-06 -1.382250e-05 -5.104706e-04
                         -5.985740e-06  3.573593e-06 -1.212431e-08  8.538922e-06
                         -1.382250e-05 -1.212431e-08  1.356678e-06  3.953835e-06
                         -5.104706e-04  8.538922e-06  3.953835e-06  1.102622e-03
                        ] rtol = 1e-6
    # Using Stata's value, R's plm has confint incorrect
    @test confint(model) ≈ [-0.0417325 0.0324457
                            -0.0072991 0.0002169
                            -0.0014674 0.0031635
                             0.0070191 0.1390408
                           ] rtol = 1e-5
    # Within
    model = fit(EconometricModel,
                @formula(CRMRTE ~ PrbConv + AvgSen + PrbPris + absorb(County)),
                data)
    #
    @test_broken nullloglikelihood(model) ≈ 2277.488 rtol = 1e-3
    @test_broken loglikelihood(model) ≈ 1643.304 rtol = 1e-3
    xtreg_fe_i = fit(EconometricModel,
                     @formula(CRMRTE ~ PrbConv + AvgSen + PrbPris + absorb(County + Year)),
                     data)
    # Random Effects
    xtreg = fit(EconometricModel,
                @formula(CRMRTE ~ PrbConv + AvgSen + PrbPris + PID(County) + TID(Year)),
                data)
    # IV Pooling
    ivreg = fit(EconometricModel,
                @formula(CRMRTE ~ PrbConv + (AvgSen ~ PrbPris)),
                data)
    # IV Between
    ivreg_be = fit(EconometricModel,
                   @formula(CRMRTE ~ PrbConv + (AvgSen ~ PrbPris) + between(County)),
                   data)
    # IV Within
    ivreg_fe = fit(EconometricModel,
                   @formula(CRMRTE ~ PrbConv + (AvgSen ~ PrbPris) + absorb(County)),
                   data)
    ivreg_fe_i = fit(EconometricModel,
                     @formula(CRMRTE ~ PrbConv + (AvgSen ~ PrbPris) + absorb(County + Year)),
                     data)
    # IV Random
    ivreg = fit(EconometricModel,
                @formula(CRMRTE ~ PrbConv + (AvgSen ~ PrbPris) + PID(County) + TID(Year)),
                data)
  end
end
# Nominal Models
@testset "Nominal Models" begin
  data = joinpath(dirname(pathof(Econometrics)), "..", "data", "insure.csv") |>
	CSV.read |>
    (data -> data[[:insure, :age, :male, :nonwhite, :site]]) |>
    (data -> dropmissing!(data, disallowmissing = true)) |>
    (data -> categorical!(data, [:insure, :site]))
  model = fit(EconometricModel,
              @formula(insure ~ 1 + age + male + nonwhite + site),
              data,
              contrasts = Dict(:insure => DummyCoding(base = "Uninsure")))
  β, V = coef(model), vcov(model)
  @test β ≈ [1.286943, 0.0077961, -0.4518496, -0.2170589, 1.211563, 0.2078123,
             1.556656, -0.0039489, 0.1098438, 0.7577178, 1.324599, -0.3801756
             ] rtol = 1e-3
  @test V ≈ Hermitian(
  			[ 0.35084518 -0.00590456 -0.02519209 -0.01837512 -0.08995794 -0.08264268  0.29928937 -0.00510107 -0.02257270 -0.01585492 -0.07392512 -0.06808314
             -0.00590456  0.00013092 -0.00046225 -0.00028066  0.00041004  0.00040470 -0.00509454  0.00011354 -0.00039387 -0.00024122  0.00033847  0.00035236
             -0.02519209 -0.00046225  0.13504644  0.01565749  0.01009559  0.00933967 -0.02288847 -0.00039497  0.11365141  0.01369872  0.00945927  0.00802593
             -0.01837512 -0.00590456 -0.02519209  0.18116609  0.01886597 -0.02821348 -0.01638629 -0.00023179  0.01314666  0.15068111  0.01664999 -0.02345786
             -0.08995794 -0.00590456 -0.02519209 -0.01837512  0.22138224 0.064050660 -0.07438947  0.00035751  0.00845363  0.01594691  0.19895500  0.05189742
             -0.08264268  0.00040470  0.00933967 -0.02821348  0.06405066 0.134170290 -0.06796097  0.00034828  0.00826890 -0.02353016  0.05182727  0.11060487
              0.29928937 -0.00509454 -0.02288847 -0.01638629 -0.07438947 -0.06796097  0.35560782 -0.00601966 -0.02683161 -0.02105510 -0.08991231 -0.08165509
             -0.00510107  0.00011354 -0.00039497 -0.00023179  0.00035751  0.00034828 -0.00601966  0.00013455 -0.00047223 -0.00028680  0.00039990  0.0004179
             -0.02257270 -0.00039387  0.11365141  0.01314666  0.00845363  0.00826890 -0.02683161 -0.00047223  0.13336253  0.01625685  0.01101870  0.00922333
             -0.01585492 -0.00024122  0.01369872  0.15068111  0.01594691 -0.02353016 -0.02105510 -0.00028680  0.01625685  0.17604391  0.02099458 -0.03008391
             -0.07392512  0.00033847  0.00945927  0.01664999  0.19895500  0.05182727 -0.08991231  0.00039990  0.01101870  0.02099458  0.22070774  0.06273376
             -0.06808314  0.00035236  0.00802593 -0.02345786  0.05189742  0.11060487 -0.08165509  0.00041790  0.00922333 -0.03008391  0.06273376  0.13899387
            ]) rtol = 1e-3
  @test sqrt.(diag(V)) ≈ [0.5923219, 0.0114418, 0.3674867, 0.4256361, 0.4705127, 0.3662926,
  						  0.5963286, 0.0115994, 0.3651883, 0.4195759, 0.4697954, 0.3728188
						 ] rtol = 1e-4
end
# Ordinal Models
@testset "Ordinal Models" begin
  data = dataset("Ecdat", "Kakadu")[[:RecParks, :Sex, :Age, :Schooling]]
  data.RecParks = convert(Vector{Int}, data.RecParks)
  data.RecParks = levels!(categorical(data.RecParks, ordered = true), collect(1:5))
  model = fit(EconometricModel,
               @formula(RecParks ~ Age + Sex + Schooling),
               data,
               contrasts = Dict(:RecParks => DummyCoding(levels = collect(1:5)))
               )
  β, V = coef(model), vcov(model)
  @test β ≈ [0.009437926, -0.015143049, -0.103911316,
             -2.92391240, -1.549266200, -0.298963900, 0.6698249] rtol = 1e-4
  @test V ≈ [ 0.00000646 -0.0000007000  0.000013250 0.00031473	0.0003161358 0.000321276 0.0003270734
             -0.00000070  0.0071633800 -0.000184930 0.00294759	0.0029647442 0.002995366 0.0029904485
              0.00001325 -0.0001849300  0.000618730 0.00288365	0.0028548209 0.002791361 0.0027318928
              0.00031473  0.0029475900  0.002883650 0.03684071	0.0289239290 0.026986237 0.026498765
              0.00031614  0.0029647400  0.002854820 0.02892393	0.0293971651 0.027160433 0.0266006675
              0.000321276 0.0029953660  0.002791361 0.02698623	0.0271604331 0.027860143 0.0270590982
              0.0003270734 0.002990448  0.002731893 0.02649876	0.0266006675 0.027059098 0.0281002038
            ] rtol = 1e-3
  @test sqrt.(diag(V)) ≈ [0.002540789, 0.084636770, 0.024874381,
                          0.191939349, 0.171456015, 0.166913580, 0.167631154] rtol = 1e-4
  data = joinpath(dirname(pathof(Econometrics)), "..", "data", "auto.csv") |>
	CSV.read |>
    (data -> data[[:rep77, :foreign, :length, :mpg]]) |>
    (data -> dropmissing!(data, disallowmissing = true))
  data.rep77 = levels!(categorical(data.rep77; ordered = true),
                       ["Poor", "Fair", "Average", "Good", "Excellent"])
  model = fit(EconometricModel,
              @formula(rep77 ~ foreign + length + mpg),
              data)
  β, V = coef(model), vcov(model)
  @test β ≈ [2.89679875, 0.08282676, 0.23076532, 17.92728, 19.86486, 22.10311, 24.69193] atol = 1e-3
  @test V ≈ [0.62578335	0.0111044038 0.012321411 2.4421745	2.4836507	2.5750058	2.6957587
             0.0111044	0.0005186744 0.001222385 0.1243671	0.1258979	0.1284727	0.1319772
             0.01232141	0.0012223847 0.004974678 0.3328169	0.3368057	0.3430257	0.3549711
             2.44217451	0.1243670884 0.332816913 30.9544666	31.0279194	31.6005528	32.5099187
             2.48365067	0.1258979463 0.336805731 31.0279194	31.4396158	31.9986912	32.9170873
             2.57500582	0.1284726666 0.343025653 31.6005528	31.9986912	32.7044311	33.6189243
             2.69575871	0.1319772302 0.354971081 32.5099187	32.9170873	33.6189243	34.8533724
            ] rtol = 1e-2
  @test sqrt.(diag(V)) ≈ [0.79106469, 0.02277442, 0.07053140,
                          5.56367384, 5.60710405, 5.71877881, 5.90367448] rtol = 1e-2
end


text:
function strip_lone_vertices(g)
    neighbors_count = map(x -> length(all_neighbors(g, x)), vertices(g))
    lone_vertex = neighbors_count .== 0
    id_mapping = (1:length(vertices(g)))[.!lone_vertex]
    new_graph, vmap = induced_subgraph(g, id_mapping)
    return new_graph
end

component_including_vertex(g, v::Union{BitArray, Array{Bool, 1}}) = map(i -> component_including_vertex(g, i), findall(v))
component_including_vertex(g, v::Int) = induced_subgraph(g, filter(x -> v in x, weakly_connected_components(g))[1])[1]


text:
"""
    get_substruct_match(mol::Mol, qmol::Mol, options::Union{Dict{String,Any},Nothing}=nothing)::Dict{String, Any}

Gets a substructure match.

```julia
smatch = get_substruct_match(mol)
```
"""
function get_substruct_match(mol::Mol, qmol::Mol, options::Union{Dict{String,Any},Nothing}=nothing)::Dict{String, Any}
    options_json::String = jsonify_details(options)
    val::Cstring = ccall((:get_substruct_match, librdkitcffi), Cstring, (Cstring, Csize_t, Cstring, Csize_t, Cstring), mol.mol[], mol.mol_size[], qmol.mol[], qmol.mol_size[], options_json)
    smatch = JSON.parse(unsafe_string_and_free(val))
    return smatch
end

"""
    get_substruct_matches(mol::Mol, qmol::Mol, options::Union{Dict{String,Any},Nothing}=nothing)::Dict{String, Any}

Gets all substructure matches.

```julia
smatches = get_substruct_matches(mol)
```
"""
function get_substruct_matches(mol::Mol, qmol::Mol, options::Union{Dict{String,Any},Nothing}=nothing)::Vector{Any}
    options_json::String = jsonify_details(options)
    val::Cstring = ccall((:get_substruct_matches, librdkitcffi), Cstring, (Cstring, Csize_t, Cstring, Csize_t, Cstring), mol.mol[], mol.mol_size[], qmol.mol[], qmol.mol_size[], options_json)
    smatches = JSON.parse(unsafe_string_and_free(val))
    return smatches
end


text:
# script to run as: julia examples/warntype.jl
using PNML, EzXML, ModuleDocstrings, InteractiveUtils

node = EzXML.ElementNode("n1/")
reg = PNML.IDRegistry()
pdict = PNML.PnmlDict()
pntd = PnmlCore()

InteractiveUtils.@code_warntype PNML.parse_file("file")
InteractiveUtils.@code_warntype PNML.parse_str("</tag>")
InteractiveUtils.@code_warntype PNML.parse_doc(EzXML.parsexml("<pnml></pnml>"))


text:
@inline function setDandΔs!(mechanism::Mechanism, diagonal::DiagonalEntry, body::Body)
    diagonal.D = ∂dyn∂vel(body, mechanism.Δt)
    diagonal.Δs = dynamics(mechanism, body)
    return
end

@inline function extendDandΔs!(mechanism::Mechanism, diagonal::DiagonalEntry, body::Body, ineqc::InequalityConstraint)
    diagonal.D += schurD(ineqc, body, mechanism.Δt) # + SMatrix{6,6,Float64,36}(1e-5*I)
    diagonal.Δs += schurf(mechanism, ineqc, body)
    return
end

@inline function setDandΔs!(mechanism::Mechanism, diagonal::DiagonalEntry{T,N}, eqc::EqualityConstraint) where {T,N}
    # diagonal.D = szeros(T, N, N)
    μ = 1e-10
    diagonal.D = -I*μ
    diagonal.Δs = g(mechanism, eqc)
    return
end

@inline function setLU!(mechanism::Mechanism, offdiagonal::OffDiagonalEntry, bodyid::Integer, eqc::EqualityConstraint)
    offdiagonal.L = ∂g∂ʳpos(mechanism, eqc, bodyid)'
    offdiagonal.U = ∂g∂ʳvel(mechanism, eqc, bodyid)
    return
end

@inline function setLU!(mechanism::Mechanism, offdiagonal::OffDiagonalEntry, eqc::EqualityConstraint, bodyid::Integer)
    offdiagonal.L = ∂g∂ʳvel(mechanism, eqc, bodyid)
    offdiagonal.U = ∂g∂ʳpos(mechanism, eqc, bodyid)'
    return
end

@inline function setLU!(offdiagonal::OffDiagonalEntry{T,N1,N2}) where {T,N1,N2}
    offdiagonal.L = szeros(T, N2, N1)
    offdiagonal.U = offdiagonal.L'
    return
end

@inline function updateLU1!(offdiagonal::OffDiagonalEntry, diagonal::DiagonalEntry, maintograndchild::OffDiagonalEntry, childtograndchild::OffDiagonalEntry)
    D = diagonal.D
    offdiagonal.L -= maintograndchild.L * D * childtograndchild.U
    offdiagonal.U -= childtograndchild.L * D * maintograndchild.U
    return
end

@inline function updateLU2!(offdiagonal::OffDiagonalEntry, diagonal::DiagonalEntry)
    Dinv = diagonal.Dinv
    offdiagonal.L = offdiagonal.L * Dinv
    offdiagonal.U = Dinv * offdiagonal.U
    return
end

@inline function updateD!(diagonal::DiagonalEntry, childdiagonal::DiagonalEntry, fillin::OffDiagonalEntry)
    diagonal.D -= fillin.L * childdiagonal.D * fillin.U
    return
end

function invertD!(diagonal::DiagonalEntry)
    diagonal.Dinv = inv(diagonal.D)
    return
end

@inline function LSol!(diagonal::DiagonalEntry, child::DiagonalEntry, fillin::OffDiagonalEntry)
    diagonal.Δs -= fillin.L * child.Δs
    return
end

function DSol!(diagonal::DiagonalEntry)
    diagonal.Δs = diagonal.Dinv * diagonal.Δs
    return
end

@inline function USol!(diagonal::DiagonalEntry, parent::DiagonalEntry, fillin::OffDiagonalEntry)
    diagonal.Δs -= fillin.U * parent.Δs
    return
end

function feasibilityStepLength!(mechanism::Mechanism)
    ldu = mechanism.ldu

    τ = 0.995
    mechanism.α = 1.

    for ineqc in mechanism.ineqconstraints
        feasibilityStepLength!(mechanism, ineqc, getineqentry(ldu, ineqc.id), τ)
    end

    return
end

function feasibilityStepLength!(mechanism, ineqc::InequalityConstraint{T,N}, ineqentry::InequalityEntry, τ) where {T,N}
    s1 = ineqc.ssol[2]
    γ1 = ineqc.γsol[2]
    Δs = ineqentry.Δs
    Δγ = ineqentry.Δγ

    for i = 1:N
        αmax = τ * s1[i] / Δs[i]
        (αmax > 0) && (αmax < mechanism.α) && (mechanism.α = αmax)
        αmax = τ * γ1[i] / -Δγ[i]
        (αmax > 0) && (αmax < mechanism.α) && (mechanism.α = αmax)
    end
    return
end

function setentries!(mechanism::Mechanism)
    graph = mechanism.graph
    ldu = mechanism.ldu

    for (id, body) in pairs(mechanism.bodies)
        isinactive(body) && continue

        for childid in directchildren(graph, id)
            setLU!(mechanism, getentry(ldu, (id, childid)), id, geteqconstraint(mechanism, childid))
        end

        diagonal = getentry(ldu, id)
        setDandΔs!(mechanism, diagonal, body)
        for childid in ineqchildren(graph, id)
            ineqc = getineqconstraint(mechanism, childid)
            calcFrictionForce!(mechanism, ineqc)
            extendDandΔs!(mechanism, diagonal, body, ineqc)
        end
    end

    for eqc in mechanism.eqconstraints
        id = eqc.id
        isinactive(eqc) && continue

        for childid in directchildren(graph, id)
            setLU!(mechanism, getentry(ldu, (id, childid)), eqc, childid)
        end

        for childid in loopchildren(graph, id)
            setLU!(getentry(ldu, (id, childid)))
        end

        diagonal = getentry(ldu, id)
        setDandΔs!(mechanism, diagonal, eqc)
    end

    return
end

function factor!(graph::Graph, ldu::SparseLDU)
    for id in graph.dfslist
        isinactive(graph, id) && continue

        sucs = successors(graph, id)
        for childid in sucs
            isinactive(graph, childid) && continue

            offdiagonal = getentry(ldu, (id, childid))
            for grandchildid in sucs
                grandchildid == childid && break
                isinactive(graph, grandchildid) && continue

                if hasdirectchild(graph, childid, grandchildid)
                    updateLU1!(offdiagonal, getentry(ldu, grandchildid), getentry(ldu, (id, grandchildid)), getentry(ldu, (childid, grandchildid)))
                end
            end
            updateLU2!(offdiagonal, getentry(ldu, childid))
        end

        diagonal = getentry(ldu, id)

        for childid in successors(graph, id)
            isinactive(graph, childid) && continue
            updateD!(diagonal, getentry(ldu, childid), getentry(ldu, (id, childid)))
        end
        invertD!(diagonal)
    end

    return
end

function solve!(mechanism::Mechanism)
    ldu = mechanism.ldu
    graph = mechanism.graph
    dfslist = graph.dfslist

    for id in dfslist
        isinactive(graph, id) && continue

        diagonal = getentry(ldu, id)

        for childid in successors(graph, id)
            isinactive(graph, childid) && continue
            LSol!(diagonal, getentry(ldu, childid), getentry(ldu, (id, childid)))
        end
    end

    for id in graph.rdfslist
        isinactive(graph, id) && continue

        diagonal = getentry(ldu, id)
        DSol!(diagonal)

        for parentid in predecessors(graph, id)
            isinactive(graph, parentid) && continue
            USol!(diagonal, getentry(ldu, parentid), getentry(ldu, (parentid, id)))
        end

        for childid in ineqchildren(graph, id)
            isinactive(graph, childid) && continue
            eliminatedsolve!(mechanism, getineqentry(ldu, childid), diagonal, id, getineqconstraint(mechanism, childid))
        end
    end

    return
end

function eliminatedsolve!(mechanism::Mechanism, ineqentry::InequalityEntry, diagonal::DiagonalEntry, bodyid::Integer, ineqc::InequalityConstraint)
    Δt = mechanism.Δt
    μ = mechanism.μ

    φ = g(mechanism, ineqc)

    Nx = ∂g∂ʳpos(mechanism, ineqc, bodyid)
    Nv = ∂g∂ʳvel(mechanism, ineqc, bodyid)

    γ1 = ineqc.γsol[2]
    s1 = ineqc.ssol[2]

    Δv = diagonal.Δs
    ineqentry.Δγ = -γ1 ./ s1 .* φ + μ ./ s1 + γ1 ./ s1 .* (Nv * Δv)
    ineqentry.Δs = s1 .- μ ./ γ1 + s1 ./ γ1 .* ineqentry.Δγ

    return
end

@inline function updatesolution!(body::Body)
    body.state.vsol[1] = body.state.vsol[2]
    body.state.ωsol[1] = body.state.ωsol[2]
    return
end

@inline function updatesolution!(eqc::EqualityConstraint)
    eqc.λsol[1] = eqc.λsol[2]
    return
end

@inline function updatesolution!(ineqc::InequalityConstraint)
    ineqc.ssol[1] = ineqc.ssol[2]
    ineqc.γsol[1] = ineqc.γsol[2]
    return
end

@inline function normΔs(body::Body)
    d1 = body.state.vsol[2] - body.state.vsol[1]
    d2 = body.state.ωsol[2] - body.state.ωsol[1]
    return dot(d1, d1) + dot(d2, d2)
end

@inline function normΔs(eqc::EqualityConstraint)
    d = eqc.λsol[2] - eqc.λsol[1]
    return dot(d, d)
end

@inline function normΔs(ineqc::InequalityConstraint)
    d1 = ineqc.ssol[2] - ineqc.ssol[1]
    d2 = ineqc.γsol[2] - ineqc.γsol[1]
    return dot(d1, d1) + dot(d2, d2)
end


text:
## %PREDICTMD_GENERATED_BY%

import PredictMDExtra
PredictMDExtra.import_all()

import PredictMD
PredictMD.import_all()

# PREDICTMD IF INCLUDE TEST STATEMENTS
import CSVFiles
import DataFrames
import FileIO
import RDatasets
import Random
import Test
# PREDICTMD ELSE
# PREDICTMD ENDIF INCLUDE TEST STATEMENTS

# PREDICTMD IF INCLUDE TEST STATEMENTS
logger = Base.CoreLogging.current_logger_for_env(Base.CoreLogging.Debug, Symbol(splitext(basename(something(@__FILE__, "nothing")))[1]), something(@__MODULE__, "nothing"))
if isnothing(logger)
    logger_stream = devnull
else
    logger_stream = logger.stream
end
# PREDICTMD ELSE
# PREDICTMD ENDIF INCLUDE TEST STATEMENTS

### Begin project-specific settings

DIRECTORY_CONTAINING_THIS_FILE = @__DIR__
PROJECT_DIRECTORY = dirname(
    joinpath(splitpath(DIRECTORY_CONTAINING_THIS_FILE)...)
    )
PROJECT_OUTPUT_DIRECTORY = joinpath(
    PROJECT_DIRECTORY,
    "output",
    )
mkpath(PROJECT_OUTPUT_DIRECTORY)
mkpath(joinpath(PROJECT_OUTPUT_DIRECTORY, "data"))
mkpath(joinpath(PROJECT_OUTPUT_DIRECTORY, "models"))
mkpath(joinpath(PROJECT_OUTPUT_DIRECTORY, "plots"))

# PREDICTMD IF INCLUDE TEST STATEMENTS
@debug("PROJECT_OUTPUT_DIRECTORY: ", PROJECT_OUTPUT_DIRECTORY,)
if PredictMD.is_travis_ci()
    PredictMD.cache_to_path!(
        ;
        from = ["cpu_examples", "breast_cancer_biopsy", "output",],
        to = [PROJECT_OUTPUT_DIRECTORY],
        )
end
# PREDICTMD ELSE
# PREDICTMD ENDIF INCLUDE TEST STATEMENTS

### End project-specific settings

### Begin logistic classifier code

Random.seed!(999)

trainingandtuning_features_df_filename = joinpath(
    PROJECT_OUTPUT_DIRECTORY,
    "data",
    "trainingandtuning_features_df.csv",
    )
trainingandtuning_labels_df_filename = joinpath(
    PROJECT_OUTPUT_DIRECTORY,
    "data",
    "trainingandtuning_labels_df.csv",
    )
testing_features_df_filename = joinpath(
    PROJECT_OUTPUT_DIRECTORY,
    "data",
    "testing_features_df.csv",
    )
testing_labels_df_filename = joinpath(
    PROJECT_OUTPUT_DIRECTORY,
    "data",
    "testing_labels_df.csv",
    )
training_features_df_filename = joinpath(
    PROJECT_OUTPUT_DIRECTORY,
    "data",
    "training_features_df.csv",
    )
training_labels_df_filename = joinpath(
    PROJECT_OUTPUT_DIRECTORY,
    "data",
    "training_labels_df.csv",
    )
tuning_features_df_filename = joinpath(
    PROJECT_OUTPUT_DIRECTORY,
    "data",
    "tuning_features_df.csv",
    )
tuning_labels_df_filename = joinpath(
    PROJECT_OUTPUT_DIRECTORY,
    "data",
    "tuning_labels_df.csv",
    )
trainingandtuning_features_df = DataFrames.DataFrame(
    FileIO.load(
        trainingandtuning_features_df_filename;
        type_detect_rows = 100,
        )
    )
trainingandtuning_labels_df = DataFrames.DataFrame(
    FileIO.load(
        trainingandtuning_labels_df_filename;
        type_detect_rows = 100,
        )
    )
testing_features_df = DataFrames.DataFrame(
    FileIO.load(
        testing_features_df_filename;
        type_detect_rows = 100,
        )
    )
testing_labels_df = DataFrames.DataFrame(
    FileIO.load(
        testing_labels_df_filename;
        type_detect_rows = 100,
        )
    )
training_features_df = DataFrames.DataFrame(
    FileIO.load(
        training_features_df_filename;
        type_detect_rows = 100,
        )
    )
training_labels_df = DataFrames.DataFrame(
    FileIO.load(
        training_labels_df_filename;
        type_detect_rows = 100,
        )
    )
tuning_features_df = DataFrames.DataFrame(
    FileIO.load(
        tuning_features_df_filename;
        type_detect_rows = 100,
        )
    )
tuning_labels_df = DataFrames.DataFrame(
    FileIO.load(
        tuning_labels_df_filename;
        type_detect_rows = 100,
        )
    )

smoted_training_features_df_filename = joinpath(
    PROJECT_OUTPUT_DIRECTORY,
    "data",
    "smoted_training_features_df.csv",
    )
smoted_training_labels_df_filename = joinpath(
    PROJECT_OUTPUT_DIRECTORY,
    "data",
    "smoted_training_labels_df.csv",
    )
smoted_training_features_df = DataFrames.DataFrame(
    FileIO.load(
        smoted_training_features_df_filename;
        type_detect_rows = 100,
        )
    )
smoted_training_labels_df = DataFrames.DataFrame(
    FileIO.load(
        smoted_training_labels_df_filename;
        type_detect_rows = 100,
        )
    )

categorical_feature_names_filename = joinpath(
    PROJECT_OUTPUT_DIRECTORY,
    "data",
    "categorical_feature_names.jld2",
    )
continuous_feature_names_filename = joinpath(
    PROJECT_OUTPUT_DIRECTORY,
    "data",
    "continuous_feature_names.jld2",
    )
categorical_feature_names = FileIO.load(
    categorical_feature_names_filename,
    "categorical_feature_names",
    )
continuous_feature_names = FileIO.load(
    continuous_feature_names_filename,
    "continuous_feature_names",
    )
feature_names = vcat(categorical_feature_names, continuous_feature_names)

single_label_name = :Class
negative_class = "benign"
positive_class = "malignant"
single_label_levels = [negative_class, positive_class]

categorical_label_names = Symbol[single_label_name]
continuous_label_names = Symbol[]
label_names = vcat(categorical_label_names, continuous_label_names)

feature_contrasts = PredictMD.generate_feature_contrasts(
    smoted_training_features_df,
    feature_names,
    )

# PREDICTMD IF INCLUDE TEST STATEMENTS
Test.@test(
    PredictMD.columns_are_linearly_independent(training_features_df)
    )
Test.@test(
    PredictMD.columns_are_linearly_independent(
        training_features_df,
        feature_names,
        )
    )
Test.@test(
    length(PredictMD.linearly_dependent_columns(df)) == 0
    )
Test.@test(
    length(
        PredictMD.linearly_dependent_columns(
            training_features_df,
            feature_names,
            ),
        ) == 0
    )
# PREDICTMD ELSE
# PREDICTMD ENDIF INCLUDE TEST STATEMENTS

show(
    logger_stream, PredictMD.linearly_dependent_columns(df)
    )

show(
    logger_stream, PredictMD.linearly_dependent_columns(
        training_features_df,
        feature_names,
        )
    )

# PREDICTMD IF INCLUDE TEST STATEMENTS
Test.@test_throws ErrorException PredictMD.singlelabelbinaryclassdataframelogisticclassifier(
    feature_names,
    single_label_name,
    single_label_levels;
    package = :thispackagedoesnotexist,
    intercept = true,
    interactions = 1,
    name = "Logistic regression",
    )
# PREDICTMD ELSE
# PREDICTMD ENDIF INCLUDE TEST STATEMENTS

logistic_classifier =
        PredictMD.singlelabelbinaryclassdataframelogisticclassifier(
        feature_names,
        single_label_name,
        single_label_levels;
        package = :GLM,
        intercept = true,
        interactions = 1,
        name = "Logistic regression",
        )

# PREDICTMD IF INCLUDE TEST STATEMENTS
PredictMD.predict(logistic_classifier, training_features_df)
PredictMD.predict_proba(logistic_classifier, training_features_df)
# PREDICTMD ELSE
# PREDICTMD ENDIF INCLUDE TEST STATEMENTS

PredictMD.fit!(logistic_classifier,
               smoted_training_features_df,
               smoted_training_labels_df) # TODO: fix this error

PredictMD.get_underlying(logistic_classifier) # TODO: fix this error

logistic_hist_training =
        PredictMD.plotsinglelabelbinaryclassifierhistogram( # TODO: fix this error
        logistic_classifier,
        smoted_training_features_df,
        smoted_training_labels_df,
        single_label_name,
        single_label_levels,
        );

# PREDICTMD IF INCLUDE TEST STATEMENTS
filename = string(
    tempname(),
    "_",
    "logistic_hist_training",
    ".pdf",
    )
rm(filename; force = true, recursive = true,)
@debug("Attempting to test that the file does not exist...", filename,)
Test.@test(!isfile(filename))
@debug("The file does not exist.", filename, isfile(filename),)
PredictMD.save_plot(filename, logistic_hist_training)
if PredictMD.is_force_test_plots()
    @debug("Attempting to test that the file exists...", filename,)
    Test.@test(isfile(filename))
    @debug("The file does exist.", filename, isfile(filename),)
end
# PREDICTMD ELSE
# PREDICTMD ENDIF INCLUDE TEST STATEMENTS

display(logistic_hist_training)
PredictMD.save_plot(
    joinpath(
        PROJECT_OUTPUT_DIRECTORY,
        "plots",
        "logistic_hist_training.pdf",
        ),
    logistic_hist_training,
    )

logistic_hist_testing =
    PredictMD.plotsinglelabelbinaryclassifierhistogram(
        logistic_classifier,
        testing_features_df,
        testing_labels_df,
        single_label_name,
        single_label_levels,
        );

# PREDICTMD IF INCLUDE TEST STATEMENTS
filename = string(
    tempname(),
    "_",
    "logistic_hist_testing",
    ".pdf",
    )
rm(filename; force = true, recursive = true,)
@debug("Attempting to test that the file does not exist...", filename,)
Test.@test(!isfile(filename))
@debug("The file does not exist.", filename, isfile(filename),)
PredictMD.save_plot(filename, logistic_hist_testing)
if PredictMD.is_force_test_plots()
    @debug("Attempting to test that the file exists...", filename,)
    Test.@test(isfile(filename))
    @debug("The file does exist.", filename, isfile(filename),)
end
# PREDICTMD ELSE
# PREDICTMD ENDIF INCLUDE TEST STATEMENTS

display(logistic_hist_testing)
PredictMD.save_plot(
    joinpath(
        PROJECT_OUTPUT_DIRECTORY,
        "plots",
        "logistic_hist_testing.pdf",
        ),
    logistic_hist_testing,
    )

show(
    logger_stream, PredictMD.singlelabelbinaryclassificationmetrics(
        logistic_classifier,
        smoted_training_features_df,
        smoted_training_labels_df,
        single_label_name,
        positive_class;
        sensitivity = 0.95,
        );
    allrows = true,
    allcols = true,
    splitcols = false,
    )

show(
    logger_stream, PredictMD.singlelabelbinaryclassificationmetrics(
        logistic_classifier,
        testing_features_df,
        testing_labels_df,
        single_label_name,
        positive_class;
        sensitivity = 0.95,
        );
    allrows = true,
    allcols = true,
    splitcols = false,
    )

logistic_calibration_curve =
    PredictMD.plot_probability_calibration_curve(
        logistic_classifier,
        smoted_training_features_df,
        smoted_training_labels_df,
        single_label_name,
        positive_class;
        window = 0.2,
        );

# PREDICTMD IF INCLUDE TEST STATEMENTS
# PREDICTMD ELSE
# PREDICTMD ENDIF INCLUDE TEST STATEMENTS

display(logistic_calibration_curve)
PredictMD.save_plot(
    joinpath(
        PROJECT_OUTPUT_DIRECTORY,
        "plots",
        "logistic_calibration_curve.pdf",
        ),
    logistic_calibration_curve,
    )

show(
    logger_stream, PredictMD.probability_calibration_metrics(
        logistic_classifier,
        testing_features_df,
        testing_labels_df,
        single_label_name,
        positive_class;
        window = 0.1,
        );
    allrows = true,
    allcols = true,
    splitcols = false,
    )

logistic_cutoffs, logistic_risk_group_prevalences =
    PredictMD.risk_score_cutoff_values(
        logistic_classifier,
        testing_features_df,
        testing_labels_df,
        single_label_name,
        positive_class;
        average_function = Statistics.mean,
        )
# PREDICTMD IF INCLUDE TEST STATEMENTS
@debug(
    string(
        "Low risk: 0 to $(logistic_cutoffs[1]).",
        " Medium risk: $(logistic_cutoffs[1]) to $(logistic_cutoffs[2]).",
        " High risk: $(logistic_cutoffs[2]) to 1.",
        )
    )
@debug(logistic_risk_group_prevalences)
# PREDICTMD ELSE
@info(
    string(
        "Low risk: 0 to $(logistic_cutoffs[1]).",
        " Medium risk: $(logistic_cutoffs[1]) to $(logistic_cutoffs[2]).",
        " High risk: $(logistic_cutoffs[2]) to 1.",
        )
    )
@info(logistic_risk_group_prevalences)
# PREDICTMD ENDIF INCLUDE TEST STATEMENTS

logistic_cutoffs, logistic_risk_group_prevalences =
    PredictMD.risk_score_cutoff_values(
        logistic_classifier,
        testing_features_df,
        testing_labels_df,
        single_label_name,
        positive_class;
        average_function = Statistics.median,
        )
# PREDICTMD IF INCLUDE TEST STATEMENTS
@debug(
    string(
        "Low risk: 0 to $(logistic_cutoffs[1]).",
        " Medium risk: $(logistic_cutoffs[1]) to $(logistic_cutoffs[2]).",
        " High risk: $(logistic_cutoffs[2]) to 1.",
        )
    )
@debug(logistic_risk_group_prevalences)
# PREDICTMD ELSE
@info(
    string(
        "Low risk: 0 to $(logistic_cutoffs[1]).",
        " Medium risk: $(logistic_cutoffs[1]) to $(logistic_cutoffs[2]).",
        " High risk: $(logistic_cutoffs[2]) to 1.",
        )
    )
@info(logistic_risk_group_prevalences)
# PREDICTMD ENDIF INCLUDE TEST STATEMENTS

logistic_classifier_filename = joinpath(
    PROJECT_OUTPUT_DIRECTORY,
    "models",
    "logistic_classifier.jld2",
    )

PredictMD.save_model(logistic_classifier_filename, logistic_classifier)

# PREDICTMD IF INCLUDE TEST STATEMENTS
logistic_classifier = nothing
Test.@test isnothing(logistic_classifier)
# PREDICTMD ELSE
# PREDICTMD ENDIF INCLUDE TEST STATEMENTS

### End logistic classifier code

# PREDICTMD IF INCLUDE TEST STATEMENTS
if PredictMD.is_travis_ci()
    PredictMD.path_to_cache!(
        ;
        to = ["cpu_examples", "breast_cancer_biopsy", "output",],
        from = [PROJECT_OUTPUT_DIRECTORY],
        )
end
# PREDICTMD ELSE
# PREDICTMD ENDIF INCLUDE TEST STATEMENTS

## %PREDICTMD_GENERATED_BY%


text:
include("mathlc.jl")
include("setup_for_centroid_julia.jl")
# @pyimport matplotlib.font_manager as fm --- marked affected locations with "--- disabled because of problem with @pyimport"

"""
Run PFVG method
# Example
```julia-repl
julia> include("PFVG.jl")
julia> runPFVG()
```
"""
function runPFVG()
    dirlc = joinpath(@__DIR__, "data/")
    lc = [dirlc*"bj.txt",dirlc*"vj.txt",dirlc*"rj.txt",dirlc*"ij.txt"]
    filter_vec = ["B","V","R","I"]
    B     = readdlm(lc[1])[:, 1:3]
    V     = readdlm(lc[2])[:, 1:3]
    R     = readdlm(lc[3])[:, 1:3]
    I     = readdlm(lc[4])[:, 1:3]

    #select same dates between B and I
    samet1,sflxB1,sflxI = mergeTimeSeries(B[:,1]',I[:,1]',B[:,2],I[:,2])
    samet1,EsflxB1,EsflxI = mergeTimeSeries(B[:,1]',I[:,1]',B[:,3],I[:,3])
    #select same dates between previous selected B and V
    samet2,sflxB2,sflxV = mergeTimeSeries(samet1',V[:,1]',sflxB1,V[:,2])
    samet2,EsflxB2,EsflxV = mergeTimeSeries(samet1',V[:,1]',EsflxB1,V[:,3])
    #select same dates between previous selected B and R
    samet3,sflxB3,sflxR = mergeTimeSeries(samet2',R[:,1]',sflxB2,R[:,2])
    samet3,EsflxB3,EsflxR = mergeTimeSeries(samet2',R[:,1]',EsflxB2,R[:,3])
    #
    samet,sflxB,sflxI   = mergeTimeSeries(samet3',samet1',sflxB2,sflxI)
    samet,EsflxB,EsflxI   = mergeTimeSeries(samet3',samet1',EsflxB2,EsflxI)

    ndata = length(samet)
    data  = zeros(4,3,ndata)

    flx   = [sflxB,sflxV,sflxR,sflxI]
    eflx  = [EsflxB,EsflxV,EsflxR,EsflxI]

    for i in 1:length(lc)
        obstime = samet #time
        flux    = flx[i]
        errflx  = eflx[i] #errflux
        data[i,1,:] = obstime
        data[i,2,:] = flux
        data[i,3,:] = errflx
    end

    #start PVFG
    # @eval using ProbabilisticFluxVariationGradient

    flx = data[:,2,:]
    sig = data[:,3,:]

    # run for one iteration only for warmup (i.e. when first run Julia pre-compiles code)
    #posterior = bpca(flx, sig,maxiter=1)
    # run for more iterations (default value of maxiter is 1000)
    #posterior = bpca(flx, sig,maxiter=500)
    posterior = bpca(flx, sig,maxiter=1000)

    galfile  = dirlc*"Galaxy_vec.txt"

    rgalfile = readdlm(galfile)[:,1]
    galaxyvec = rgalfile[:,1]

    g     = galaxyvec #galaxy fluxes in mJy
    randx = noisyintersectionvi(posterior = posterior, g = g)

    #return random galaxy values
    nsamp = 10000
    rangal = zeros(4,nsamp) #(filters,column,numbers)

    for i in 1:length(lc)
        for j in 1:nsamp
            galrandom = randx()[i]
            rangal[i,j] = galrandom
        end
        #save to file
        writedlm(dirlc*"PFVG.dist.3C120."*filter_vec[i]*".txt", rangal[i,:])
    end
end

"""
Plot PFVG host-galaxy distributions
# Example
```julia-repl
julia> include("PFVG.jl")
julia> PlotPFVG3C120()
```
"""
function PlotPFVGdist()
    dirlc = joinpath(@__DIR__, "data/")
    distfiles = [dirlc*"PFVG.dist.3C120.B.txt",dirlc*"PFVG.dist.3C120.V.txt",dirlc*"PFVG.dist.3C120.R.txt",dirlc*"PFVG.dist.3C120.I.txt"]
    color_vec  = ["blue","green","brown","red"]
    filter_vec = ["B","V","R","I"]

    fvggal = [1.69;3.89;5.90;8.97] # FVG results for 3C120
    snvec  = [46;82;87;112] #Light curves S/N

    # prop = fm.FontProperties(size=12) --- disabled because of problem with @pyimport
    pygui(true)
    fig  = figure(1,figsize=(15,15))

    for i in 1:length(distfiles)
        distmp   = readdlm(distfiles[i])[:, 1]
        distfvg  = distmp[:,1]
        tmp = 220+i
        ax = subplot(tmp)
        binsize=50
        centgal,centgal_loerr,centgal_uperr = centerror(distfvg)
        text1 = L"$f^{host}$ = "*@sprintf("%.3f", centgal)*"+"*@sprintf("%.3f", centgal_loerr)*"-"*@sprintf("%.3f", centgal_uperr)
        ax.hist(distfvg,bins=binsize,alpha=0.8,color=color_vec[i],label=text1)
        ax.fill_between([centgal-centgal_loerr, centgal+centgal_uperr],[0, 0],[1.1*150, 1.1*150], linestyle="--", color="black", alpha=1.0)
        ax.plot([fvggal[i],fvggal[i]],[0,1000],"--",linewidth=2,color="black")
        ax.text(centgal-0.008*centgal,800,"S/N = "*@sprintf("%.1f", snvec[i]),color="black",fontsize=15)

        if i==1
            ax.text(0.9,800,"3C120",fontsize=15)
            #ax.set_xlim(0.6,1.9)
        end
        ax.set_ylim(0,1000)
        ax.set_xlabel(L"$f_{host}$ "*filter_vec[i]*" [mJy]",fontsize=15)
        ax.set_ylabel("Frequency",fontsize=15)
        ax.tick_params(direction="in",labelsize=14)
        ax.legend(loc="upper left")#,prop=prop) --- disabled because of problem with @pyimport
    end
    savefig(dirlc*"PFVG3C120.pdf",bbox_inches="tight")
end

"""
Inspect 3C120 light curves
# Example
```julia-repl
julia> include("PFVG.jl")
julia> Check3C120()
```
"""

function Check3C120()
    dirlc = joinpath(@__DIR__, "data/")
    #read and plot lc's
    lc = [dirlc*"bj.txt",dirlc*"vj.txt",dirlc*"rj.txt",dirlc*"ij.txt"]
    filtervec = ["B","V","R","I"]
    colorvec  = ["blue","green","brown","red"]

    #find the minimum day for all the observations:
    minmjd = zeros(length(lc))
    for i in 1:length(lc)
        glc    = readdlm(lc[i])[:, 1:3]
        mjd    = glc[:,1] #time in MJD
        minmjd[i] = minimum(mjd)
    end
    mmjd = minimum(minmjd)

    #read separate
    B     = readdlm(lc[1])[:, 1:3]
    V     = readdlm(lc[2])[:, 1:3]
    R     = readdlm(lc[3])[:, 1:3]
    I     = readdlm(lc[4])[:, 1:3]

    #select same dates between B and I
    samet1,sflxB1,sflxI = mergeTimeSeries(B[:,1]',I[:,1]',B[:,2],I[:,2])
    samet1,EsflxB1,EsflxI = mergeTimeSeries(B[:,1]',I[:,1]',B[:,3],I[:,3])
    #select same dates between previous selected B and V
    samet2,sflxB2,sflxV = mergeTimeSeries(samet1',V[:,1]',sflxB1,V[:,2])
    samet2,EsflxB2,EsflxV = mergeTimeSeries(samet1',V[:,1]',EsflxB1,V[:,3])
    #select same dates between previous selected B and R
    samet3,sflxB3,sflxR = mergeTimeSeries(samet2',R[:,1]',sflxB2,R[:,2])
    samet3,EsflxB3,EsflxR = mergeTimeSeries(samet2',R[:,1]',EsflxB2,R[:,3])
    #
    samet,sflxB,sflxI   = mergeTimeSeries(samet3',samet1',sflxB2,sflxI)
    samet,EsflxB,EsflxI   = mergeTimeSeries(samet3',samet1',EsflxB2,EsflxI)

    samet= samet .- mmjd

    #ploting
    # prop = fm.FontProperties(size=15)
    figure(1)
    clf()

    for i in 1:length(lc)
        glc     = readdlm(lc[i])[:, 1:3]
        mjd     = glc[:,1] #time
        flux    = glc[:,2] #flux
        errflx  = glc[:,3] #errflux
        obsdays = mjd .- mmjd

        plot(obsdays,flux,ls="", marker=".", ms=4,color=colorvec[i],label=filtervec[i])
        errorbar(obsdays, flux, yerr=errflx, fmt="o", ms=4, color=colorvec[i])
        legend(loc="upper right")#,prop=prop) --- disabled because of problem with @pyimport

    end
    plot(samet,sflxB,ls="", marker="o", ms=8,color="black")
    plot(samet,sflxV,ls="", marker="o", ms=8,color="black")
    plot(samet,sflxR,ls="", marker="o", ms=8,color="black")
    plot(samet,sflxI,ls="", marker="o", ms=8,color="black")
    #xlim(0.0,maxday)
    #ylim(3.0,8.0)
    xlabel("Days after MJD"*string(mmjd),fontsize=20)
    ylabel("Flux (mJy)",fontsize=20)
    tick_params(direction="in",labelsize=14)
    gcf()
end


text:
exists(lib::Ptr, sym::Symbol) = Libdl.dlsym_e(lib, sym) != C_NULL

function udp_srv(port::Int)
    r = Future()

    sock = UDPSocket()
    bind(sock, ip"127.0.0.1", port)

    @async begin
        put!(r, String(recv(sock)))
        close(sock)
    end

    return r
end

function tcp_srv(port::Int)
    r = Future()

    server = listen(ip"127.0.0.1", port)
    @async begin
        while !isready(r)
            sock = accept(server)
            put!(r, readuntil(sock, '\0'; keep=true))
        end
        close(server)
    end

    return r
end


text:
import PyCall: PyObject

const XFloat  = Union{Float64,Float32,Float16}
const XInt    = Union{Int64,Int32,Int16}
const XScalar = Union{XFloat,XInt}

const julia2numpy = Dict(
  Float64 => "float64",
  Float32 => "float32",
  Int64   => "int64",
  Bool    => "bool")

const numpy2julia = Dict(v => k for (k, v) in julia2numpy)

numpytype(T) = xlaclient.np.dtype(julia2numpy[T])
primitivetype(T) = xlaclient.dtype_to_etype(numpytype(T))
juliatype(T) = numpy2julia[T.name]

# Shapes

struct XShape{T,N}
  dims::NTuple{N,Int}
end

XShape(T::Type{<:Union{XScalar,Bool}}, sh::NTuple{N,Integer}) where N = XShape{T,N}(sh)

shapeof(x::AbstractArray) = XShape{eltype(x),ndims(x)}(size(x))

Base.eltype(::XShape{T}) where T = T
Base.ndims(::XShape{T,N}) where {T,N} = N

PyObject(sh::XShape) = xlaclient.Shape.array_shape(numpytype(eltype(sh)), sh.dims)
pyshape(sh::XShape) = PyObject(sh)

Base.show(io::IO, sh::XShape) = print(io, eltype(sh), "[", join(sh.dims, ","), "]")

function XShape(sh::PyObject)
  T = juliatype(sh.numpy_dtype())
  size = sh.dimensions()
  XShape{T,length(size)}(size)
end

shapeof(p::PyObject) = p.is_array() ? XShape(p) : (shapeof.(p.tuple_shapes())...,)

pyshape(x::Tuple) = xlaclient.Shape.tuple_shape(pyshape.(x))

pyshape(x::Type{<:Union{XScalar,Bool}}) = pyshape(XShape(x, ()))

# Values

struct XArray{T,N} <: AbstractArray{T,N}
  buffer::PyObject
end

function setup_finaliser(x)
  delete = x.buffer.delete # work around a segfault on exit
  finalizer(buf -> ispynull(delete) || delete(), x.buffer)
end

function XArray(data::Array{<:XScalar})
  buffer = xlaclient.Buffer.from_pyval(data)
  x = XArray{eltype(data),ndims(data)}(buffer)
  setup_finaliser(x)
  return x
end

function XArray(buf::PyObject, own = true)
  sh = XShape(buf.shape())
  x = XArray{eltype(sh),ndims(sh)}(buf)
  own && setup_finaliser(x)
  return x
end

PyObject(x::XArray) = x.buffer

# PyCall seems to sometimes convert arrays, sometimes not.
npy2julia(x::PyObject) = reshape(x.transpose().flatten(), x.shape)
npy2julia(x) = x

Base.size(x::XArray) = x.buffer.shape().dimensions()
Base.collect(x::XArray) = npy2julia(x.buffer.to_py())
Base.print_array(io::IO, x::XArray) = Base.print_array(io, collect(x))
Base.show_vector(io::IO, x::XArray) = Base.show_vector(io, collect(x))
Base._show_nonempty(io::IO, x::XArray, prefix::String) = Base._show_nonempty(io, collect(x), prefix)

scalar(x::PyObject) = get(x, ())
scalar(x::Array) = x[]
scalar(x::XArray{T,0}) where T = scalar(x.buffer.to_py()) # Array or PyObject? Seems to be random
scalar(x::XArray) = x

xla(x::XArray) = x
xla(x::AbstractArray{<:XScalar}) = XArray(x)
xla(x::Number) = XArray(fill(x))
xla(x::Tuple) = xla.(x)

function wrapvalue(p::PyObject)
  p.shape().is_tuple() ? (wrapvalue.(p.destructure())...,) : scalar(XArray(p))
end

default_device() = xlaclient.get_local_backend().devices()[1]

buffer(x::Array{<:XScalar}) = xlaclient.Buffer.from_pyval(x)
buffer(x::XScalar) = xlaclient.Buffer.from_pyval(x)
buffer(x::XArray) = x.buffer

# IR Builder

shapeof(builder, op) = shapeof(builder.GetShape(op))

function settypes!(builder, comp::IR, ops...; with = identity)
  Ts = map(op -> with(shapeof(builder, op)), ops)
  argtypes(comp)[:] = [Ts...]
  return comp
end

function build(ir::IR)
  builder = xlaclient.XlaBuilder("")
  env = Dict()
  resolve(x::Variable) = env[x]
  resolve(x::QuoteNode) = const!(builder, x.value)
  resolve(x) = const!(builder, x)
  for (v, T, i) in zip(arguments(ir), argtypes(ir), 1:length(arguments(ir)))
    env[v] = xlaclient.ops.Parameter(builder, i-1, pyshape(T))
  end
  for (v, st) in ir
    ex = st.expr
    if isexpr(ex, :call)
      env[v] = build!(builder, ex.args[1], resolve.(ex.args[2:end])...)
    elseif ex isa IR
      env[v] = ex
    elseif isexpr(ex)
      error("Invalid XLA expression $(ex)")
    elseif ex isa Variable
      env[v] = env[ex]
    else
      env[v] = const!(builder, ex)
    end
  end
  if isreturn(blocks(ir)[end])
    ret = returnvalue(blocks(ir)[end])
    # TODO handle the variable case
    !(ret isa Variable) && const!(builder, ret)
  end
  return builder.Build()
end

function compile(ir::IR)
  ir = controlflow(ir)
  comp = xlaclient.get_local_backend().compile(build(ir))
  return (xs...) -> wrapvalue.(comp.Execute(buffer.(xs)))
end


text:
# ## 多項式回帰
# ### 例52
# sin関数に乱数を足した観測データを多項式で回帰します。
# \lineskip
# データの生成
using Plots, Joe, Random
using Plots.PlotMeasures # hide
Plots.reset_defaults() # hide
gr() # hide
default( # hide
    titlefont  = font("JuliaMono", default(:titlefontsize),  ), # hide
    guidefont  = font("JuliaMono",  default(:guidefontsize),  ), # hide
    tickfont   = font("JuliaMono", default(:tickfontsize),   ), # hide
    legendfont = font("JuliaMono",  default(:legendfontsize), ), # hide
    left_margin = 30px, # hide
    bottom_margin = 30px # hide
) # hide
Random.seed!(12)
n = 100; x = randn(n); y = sin.(x) + randn(n);
p52 = scatter(x,y,xlabel="x", ylabel="y", label=false, legend=:topleft);
# 3、5、7次の多項式で回帰します。
# 多項式回帰用にデータセットを作る関数を定義したので、これを使ってみます。
# ```julia
# polynomial(x::Vector, P::Int) = hcat([x.^p for p in 1:P]...)
# polynomial(x::AbstractRange, P::Int)= polynomial(collect(x),P)
# ```
using Joe:polynomial
p_set = [3,5,7]
x_seq = -3:0.1:3
for p in p_set
    global x, y, x_seq, p52 # hide
    X = polynomial(x,p)
    β̂ = multiple_regression(X,y)
    ŷ = insert_ones(polynomial(x_seq,p))*β̂
    plot!(p52, x_seq,ŷ,label="p = $p")
end
savefig(p52,joinpath(@OUTPUT,"fig6-1.svg")) # hide
# \fig{fig6-1}
# ### 例53
# 矩形波に雑音をノイズをのせたデータを作って、フィッティングします。
# まずはデータを生成します。
using Plots, Random
Random.seed!(123)
n = 100; x = randn(n) * π
y = abs.(round.(x) .%2) * 2  .- 1+ randn(n)*0.2
p53 = scatter(x,y,label= false);
# 基底関数が、$f_1(x) = 1$、 $f_2(x) = \cos(x) $、 $f_3(x) = \cos(2x)$、
# $f_4(x) = \cos(3x)$の場合と、$f_1(x) = 1$、 $f_2(x) = \sin(x) $、 $f_3(x) = \sin(2x)$、
# $f_4(x) = \sin(4x)$の場合とでそれぞれ回帰分析します。
# 準備として先ほどのpolynomial関数をべき乗以外にも対応できるように多重ディスパッチで拡張します。
# ```julia
# polynomial(x::Vector, P::Int, f::Function) = hcat([f.(p*x) for p in 1:P]...)
# polynomial(x::AbstractRange, P::Int, f::Function) = polynomial(collect(x), P, f)
# ```
# $(-8,8)$の区間でフィッティングします。
x_seq = -8:0.2:8
# cosの場合
X_cos = polynomial(x,3,cos)
β̂_cos = multiple_regression(X_cos,y)
ŷ_cos = insert_ones(polynomial(x_seq,3,cos))*β̂_cos;
# sinの場合
X_sin = polynomial(x,3,sin)
β̂_sin = multiple_regression(X_sin,y)
ŷ_sin = insert_ones(polynomial(x_seq,3,sin))*β̂_sin;
# 可視化します。
plot!(p53,x_seq,ŷ_cos,label="cos")
plot!(p53,x_seq,ŷ_sin,label="sin")
savefig(p53,joinpath(@OUTPUT,"fig6-2.svg")) # hide
# \fig{fig6-2}
# \lineskip
# cosを基底にとったほうが、上手くフィッティングできることが確認できました。
# ## スプライン回帰
# スプライン回帰については、まずスプラインに関する複合型を作って、
# function-like objectを実装しました。
# spline_matrix関数は、あとで変更できるようにキーワード引数にしました。
#
# ```julia
# using Parameters
# @with_kw mutable struct Spline{T<:Number}
#     xmin::T = -1.0
#     xmax::T = 1.0
#     K::Int = 5
# end
# function spline_matrix(x::Vector,K::Int,xmin,xmax)
#     n = length(x)
#     Knots = range(xmin,stop=xmax,length=K)
#     X = zeros(n,K+4)
#     X[:,1] .= 1
#     X[:,2] .= x
#     X[:,3] .= x.^2
#     X[:,4] .= x.^3
#     for j in 1:K
#         X[:,j+4] = @. max(x - Knots[j], 0)^3
#     end
#     return X
# end
# function (s::Spline)(x::Vector, y::Vector, x_pred::Vector;
#           spline_matrix::Function=spline_matrix)
#     @unpack xmin, xmax, K = s
#     Knots = range(xmin,stop=xmax,length=K)
#     X = spline_matrix(x,K,xmin,xmax)
#     X_pred = spline_matrix(x_pred,K,xmin,xmax)
#     β = (X'X)\X'y
#     ypred = X_pred*β
# end
# (s::Spline)(x::Vector, y::Vector, x_pred::AbstractRange;kwargs...) = s(x,y,collect(x_pred);kwargs...)
# ```
# ### 例54
# sinカーブに乱数を加えたデータにスプライン回帰でフィッティングします。
# 分割点数が7以上でデータに追随していることがわかります。
using Plots, Random
using Joe: Spline
using Joe
Random.seed!(123)
n = 100; x = randn(n)*2π; y = sin.(x) + 0.2randn(n)
p54 = scatter(x,y,xlims=(-5,5),xlabel="x",ylabel="f(x)",label=false)
K_set = 5:2:9
for K in K_set
    spline = Spline(xmin=-2π, xmax=2π, K=K)
    u_seq = -5:0.2:5
    v_seq = spline(x,y,u_seq)
    plot!(p54, u_seq,v_seq, label="K = $K")
end
savefig(p54,joinpath(@OUTPUT,"fig6-5.svg")) # hide
# \fig{fig6-5}
# ### 例56 自然なスプライン曲線
# スプライン曲線の両端が直線になっているものです。
# 原著では、関数dとか関数hを定義して使っていました。
# そんなことはしたくないけど、いい名前もないので、関数内で定義して使うことにしました。
# ```julia
# function natural_spline_matrix(x::Vector, K::Int, xmin, xmax)
#     n = length(x)
#     Knots = range(xmin,stop=xmax,length=K)
#     X = zeros(n,K)
#     X[:,1] .= 1
#     X[:,2] .= x
#     d(x::Vector,a::Number,b::Number) =  @. (max(x-a,0)^3 - max(x-b,0)^3) / (b-a)
#     for j in 1:K-2
#         X[:,j+2] = d(x,Knots[j],Knots[end]) - d(x,Knots[end-1],Knots[end])
#     end
#     return X
# end
# ```
# 例54と同じデータに、ただのスプラインと自然なスプラインの回帰を適用し、
# 結果を可視化します。まずデータの生成とスプライン型を定義します。
using Plots, Random
using Joe: Spline
using Joe
Random.seed!(123)
n = 100; x = randn(n)*2π; y = sin.(x) + 0.2randn(n)
spline6 = Spline(xmin=-5, xmax=5, K=6);
spline11 = Spline(xmin=-5, xmax=5, K=11);
# 次に回帰を行います。自然なスプラインで回帰を行う場合には、キーワード引数を指定します。
u_seq = -6:0.02:6;
v_seq6 = spline6(x,y,u_seq);
v_seq11 = spline11(x,y,u_seq);
w_seq6 = spline6(x,y,u_seq;splinematrix=Joe.natural_spline_matrix);
w_seq11 = spline11(x,y,u_seq;splinematrix=Joe.natural_spline_matrix);
# 可視化します。
p56_1 = scatter(x,y,xlims=(-7,7),xlabel="x",ylabel="f(x),g(x)",label=false,title="K=6")
p56_2 = scatter(x,y,xlims=(-7,7),xlabel="x",ylabel="f(x),g(x)",label=false,title="K=11")
plot!(p56_1,u_seq,v_seq6,label="spline")
plot!(p56_1,u_seq,w_seq6,label="natural spline")
vline!(p56_1,[-5,5];linewidth=1,label=false)
vline!(p56_1,range(-5,stop=5,length=6);linestyle=:dash,label=false)
plot!(p56_2,u_seq,v_seq11,label="spline")
plot!(p56_2,u_seq,w_seq11,label="natural spline")
vline!(p56_2,[-5,5];linewidth=1,label=false)
vline!(p56_2,range(-5,stop=5,length=11);linestyle=:dash,label=false)
p56 = plot(p56_1,p56_2,layout=(2,1))
savefig(p56,joinpath(@OUTPUT,"fig6-7.svg")) # hide
# \fig{fig6-7}
# ## 平滑化スプライン
# 平滑化スプラインは区切り点を指定するときにデータ点自身を使っていたので、
# 以下のように書き直しました。(function-like objectにするまでもなかった。)
# ```julia
# struct SmoothingSpline
#     x::Vector
# end
# function (s::SmoothingSpline)(y::Vector,x_pred::Union{Vector,AbstractRange};λ=0)
#     x = s.x
#     X = smoothing_spline_matrix(x,x)
#     G = green_silverman(x)
#     γ = (X'X + λ*G)\X'y
#     X_pred = smoothing_spline_matrix(x_pred,x)
#     ypred = X_pred*γ
# end
# ```
# ただし、smoothing_spline_matrixとgreen_silverman関数は以下の通りです。
# ```julia
# function smoothing_spline_matrix(x::Vector,Knots::Union{Vector, AbstractRange})
#     n = length(x)
#     K = length(Vector)
#     X = zeros(n,K)
#     X[:,1] .= 1
#     X[:,2] .= x
#     d(x::Vector,a::Number,b::Number) =  @. (max(x-a,0)^3 - max(x-b,0)^3) / (b-a)
#     for j in 1:K-2
#         X[:,j+2] = d(x,Knots[j],Knots[end]) - d(x,Knots[end-1],Knots[end])
#     end
#     return X
# end

# function green_silverman(x::Vector)
#     n = length(x)
#     g = zeros(n,n)
#     for i in 3:n ,j in i:n
#         g[j,i] =(12(x[end]-x[end-1])*(x[end-1]-x[j-2])*(x[end-1]-x[i-2])+
#                 (x[end-1]-x[j-2])^2*(12x[end-1]+6x[j-2]-18x[i-2]))/
#                 (x[end]-x[i-2])/(x[end]-x[j-2])
#         g[i,j] = g[j,i]
#     end
#     return g
# end
# ```
# ### 例57
# $x\in(-5,5)$の一様分布からサンプルを生成します。
# 原著では、サンプル生成コードが間違っているようなので、ipynbの記述を参照しました。
using Plots, Random, Distributions
using Joe: SmoothingSpline
using Joe
Random.seed!(11)
n = 100;
x = rand(Uniform(-5,5),n); y = x .+ 2sin.(x) + randn(n)
index = sortperm(x); x = x[index]; y = y[index]

p57 = scatter(x,y,label=false)

sspline = SmoothingSpline(x)
u_seq = -8:0.02:8
λ_set = [40,400,1000]
for λ in λ_set
    global u_seq, p57 # hide
    v_seq = sspline(y,u_seq,λ=λ)
    plot!(p57,u_seq,v_seq,label="λ = $λ")
end
p57
plot!(p57,u_seq,sspline(y,u_seq,λ=1),label="λ = 1")
savefig(p57,joinpath(@OUTPUT, "fig6-8.svg")) # hide
# \fig{fig6-8}
# \lineskip
# $\lambda$が大きくなると、直線に近くなっていくようです。
# ### 例58 クロスバリデーションによる平滑化スプラインの正則化乗数の最適化
# 3章で使った線形回帰の場合の高速なクロスバリデーションの公式
# $$CV[\lambda] \coloneqq \sum_{S}||(I - H_S[\lambda])^{-1}e_S ||^2$$
# $$H_S[\lambda] \coloneqq X_s(X^\mathsf{T}X+\lambda G)^{-1}X^\mathsf{T}_S$$
# を利用して、$\lambda$の最適化を行います。
# 以下のcv_ss_fast関数を定義して利用します。
# ```julia
# function cv_ss_fast(X::AbstractArray,y::Vector,G::AbstractArray,K::Int;λ=0)
#     n = length(y); m = round(Int,n/K)
#     H = X*((X'X + λ*G)\X')
#     e = (I(n) - H)*y
#     S = 0
#     @views for j in 1:K
#         test = j*m-m+1:j*m;
#         err = (I(m)-H[test,test]) \ e[test]
#         S += dot(err,err)
#     end
#     return S/n,tr(H)
# end
# ```
# 3章の使ったcv_fast関数とほぼ同様ですが、定数項の列を付けないことと、
# 正則化乗数が引数になっていること、クロスバリデーション誤差だけでなく、
# Hのトレースも返していることなどが異なります。
# まずはデータを生成します。
using Random, Distributions
Random.seed!(11)
n = 100; x = rand(Uniform(-5,5),n) ; y = x -0.02sin.(x)-0.1randn(n);
index = sortperm(x); x = x[index]; y = y[index];
# XとGを計算します。
X = Joe.smoothing_spline_matrix(x,x);
G = Joe.green_silverman(x);
# クロスバリデーション誤差とHのトレース(有効自由度)を求めます。
result = hcat([collect(Joe.cv_ss_fast(X,y,G,n,λ=λ)) for λ in 1:50]...)
using Plots
p58 = plot(result[2,:], result[1,:],label= false,
            xlabel="有効自由度", ylabel="CVによる予測誤差",
            title = "有効自由度とCVによる予測誤差")
savefig(p58,joinpath(@OUTPUT,"fig6-9.svg"))
# \fig{fig6-9}
# \lineskip
# (注)ちなみに生成するデータによって、グラフの形はかなり変わります。
# ## 局所回帰
# ### 例60 nadaraya-watson推定量
# nadaraya-watson推定量の計算の実装です。
# カーネルは関数の引数にして、別のものに変更できるようにしました。(ここでは他のカーネルは出てきませんが。)
# ```julia
# epanechnikov(x,y,λ) = norm(x-y)/λ |> c -> max(0.75*(1-c^2),0)
#
# function nadaraya_watson_estimator(x_observed::Vector,y_observed::Vector,kk::Function,x_pred::TP;λ=1.0) where {TP<:Number}
#     n = length(y_observed)
#     S = 0; T = 0;
#     for i in 1:n
#         S += kk(x_observed[i],x_pred,λ)*y_observed[i]
#         T += kk(x_observed[i],x_pred,λ)
#     end
#     result = T == 0 ? 0 : S/T
#     return result
# end
# ベクトル化にも対応しています。
# function nadaraya_watson_estimator(x_observed::Vector,y_observed::Vector,kk::Function,x_pred::Union{Vector,AbstractRange};λ=1.0)
#     y_pred = [nadaraya_watson_estimator(x_observed,y_observed,kk,xp; λ) for xp in x_pred]
# end
# ```
# \note{epanechnikoc関数に絶対値の処理は必要になりますが、abs関数を使うよりもLinearAlgebra.jlのnormを使った方が、多次元化に対応できて良いと思います。}
# トイデータとしてsin関数にノイズをのせたデータを作成します。
using Random
Random.seed!(123)
n=250;x = 2randn(n); y = sin.(2π*x) + randn(n)/4

using Plots
p58 = scatter(x,y,xlims=(-3,3),label=false)
xx = -3:0.1:3
using Joe:nadaraya_watson_estimator, epanechnikov
using LinearAlgebra
# $\lambda=0.05,0.25$の場合にプロットしてみます。
plot!(p58,xx,nadaraya_watson_estimator(x,y,epanechnikov,xx;λ=0.05),
         label="λ = 0.05");
plot!(p58,xx,nadaraya_watson_estimator(x,y,epanechnikov,xx;λ=0.25),
         label="λ = 0.25");
# クロスバリデーションで$\lambda$を最適化してプロットしてみます。
λ_seq = 0.05:0.01:1
SS_min = Inf
λ_best = λ_seq[1]
for λ in λ_seq
    global SS_min, λ_best,n,x,y # hide
    SS = 0
    m = Int(n/10)
    for k in 1:10
        test = k*m-m+1:k*m
        train = setdiff(1:n,test)
        y_pred = nadaraya_watson_estimator(x[train],y[train],epanechnikov,x[test];λ=λ)
        SS += dot(y[test]-y_pred,y[test]-y_pred)
    end
    if SS < SS_min
        SS_min = SS
        λ_best = λ
    end
end

plot!(p58,xx,nadaraya_watson_estimator(x,y,epanechnikov,xx;λ=λ_best), label="λ = λ_best")
savefig(p58,joinpath(@OUTPUT,"fig6-10.svg")) # hide
# \fig{fig6-10}
# \lineskip
# $\lambda_{best}$はこの程度の値になるようです。
@show λ_best;
# ### 例61 局所線形回帰
# カーネル関数$k$を成分に持った対角行列$\bm{W}$を
# $$\bm{W} = \begin{bmatrix}
#  k(x,x_1) & \cdots & 0 \\
#  \vdots & \ddots & \vdots \\
#  0 & \cdots & k(x,x_n) \end{bmatrix} $$
# とした時、
# $$ L \coloneqq (\bm{y} - \bm{Xβ}(x))^\mathsf{T}\bm{W}(\bm{y} - \bm{Xβ}(x)) $$
# を最小化する$\bm{\beta}(x)$を求めます。
# $$ \dfrac{\partial L}{\partial \bm{\beta}} =
# -2\bm{X}^\mathsf{T}\bm{W}(\bm{y} - \bm{Xβ}(x))=0$$とおくと、
# $$ \bm{\beta}(x) = (\bm{X}^\mathsf{T}\bm{WX})^{-1}\bm{X}^\mathsf{T}\bm{Wy}$$
# で係数を求めることができます。この処理を関数で実装するわけですが、原著のlocalという関数名は良くないので、
# local_regressionという名前で実装しました。例では1次元に限定していましたが、
# 多次元データにも対応させています。
# ```julia
# function local_regression(x::Matrix{T},
#                           y::Vector{T},
#                           kernel::Function;
#                           x_pred=x, λ=1) where {T<:Number}
#     N,_ = size(x)
#     @assert N==length(y)
#     X = insert_ones(x)
#     m,_ = size(x_pred)
#     y_pred = Vector{T}(undef,m)
#     for j in 1:m
#         W = [kernel(x[i,:],x_pred[j,:],λ) for i in 1:N] |> diagm
#         β̂ = (X'*W*X)\X'*W*y
#         y_pred[j] = dot(insert_ones(x_pred[j,:]),β̂)
#     end
#     return y_pred
# end
# function local_regression(x::Vector{T}, y::Vector{T},kernel::Function;x_pred=x,λ=1) where {T<:Number}
#     return local_regression(x[:,:],y,kernel;x_pred=x_pred[:,:],λ=λ)
# end
# ```
# トイデータとしてsin関数にノイズをのせたデータを作成します。
using Random,Joe
Random.seed!(123)
n=30; x = 2π*rand(n) .- π; y = sin.(x) .+ randn()
p61 = scatter(x,y,label=false)
m = 200;U = -π:π/m:π;
# 正則化パラメータ$\lambda=1$の条件で回帰分析を行ってプロットします。
V = Joe.local_regression(x,y,Joe.epanechnikov;x_pred=U)
plot!(p61,U,V,label=false,title = "局所線形回帰(p=1, N=30)")
savefig(p61,joinpath(@OUTPUT,"fig6-11.svg")) # hide
# \fig{fig6-11}
# \lineskip
# 予測したい場所の数だけ回帰分析の計算を行わないといけないので、元データの数が多い時は
# 計算が大変そうです。
# ## 一般化加法モデル
# ### 例62
# 基底関数の数が有限個なら線形回帰の手法で係数を求められるけど、
# 数が多すぎると逆行列を解くのが大変でどうしましょうとい３う話。
# バックフィッティングは例63で実装します。
# ### 例63
# 以下のようなpolyfit関数を実装しました。
# ```julia
# function polyfit(x,y;x_pred=x, P::Int=3)
#     X = polynomial(x,P)
#     X_pred=polynomial(x_pred,P)
#     β̂ = multiple_regression(X,y)
#     y_pred = insert_ones(X_pred) * β̂
# end
# ```
# [例６１](#例61)のデータに一般加法モデルを適用します。
# 先に多項式でフィッティングして残りを局所線形回帰で求めていきます。
using Random,Joe
Random.seed!(123)
n=30; x = 2π*rand(n) .- π; sort!(x); y = sin.(x) .+ randn()
y₁ = zeros(n); y₂ = zeros(n)
for k in 1:10
    global x,y,y₁,y₂ # hide
    y₁ = Joe.polyfit(x,y-y₂)
    y₂ = Joe.local_regression(x,y-y₁,Joe.epanechnikov)
end
using Plots
p62_1 = plot(x,y₁,label=false,xlabel="x",ylabel="f(x)",title="多項式回帰(3次)")
p62_2 = plot(x,y₂,label=false,xlabel="x",ylabel="f(x)",title="局所線形回帰")
p62　= plot(p62_1,p62_2,layout=(1,2))
savefig(p62,joinpath(@OUTPUT,"fig6-12.svg")) # hide
# \fig{fig6-12}
# \lineskip
# 大体３次多項式で回帰されていて、残りが局所線形回帰されていることが分かります。
# 先に3次多項式で回帰したためでしょうか。 確認してみることにします。
using Random,Joe
Random.seed!(123)
n=30; x = 2π*rand(n) .- π; sort!(x); y = sin.(x) .+ randn()
y₁ = zeros(n); y₂ = zeros(n)
using Plots
for k in 1:10
    global x,y,y₁,y₂ # hide
    y₂ = Joe.local_regression(x,y-y₁,Joe.epanechnikov)
    y₁ = Joe.polyfit(x,y-y₂)
end
p62_3 = plot(x,y₁,label=false,xlabel="x",ylabel="f(x)",title="多項式回帰(3次)")
p62_4 = plot(x,y₂,label=false,xlabel="x",ylabel="f(x)",title="局所線形回帰")
p62_5　= plot(p62_3,p62_4,layout=(1,2))
savefig(p62_5,joinpath(@OUTPUT,"fig6-12-2.svg")) # hide
# \fig{fig6-12-2}
# \lineskip
# 多項式回帰の影響が減って、局所線形回帰の成分が大きくなることが確認できました。


text:
using Test
using SDPT3

# Example page 4 of SeDuMi_Guide_105R5.pdf
@testset "Linear Programming example" begin
    blk = ["l"  4]
    c = [ 1.0
         -1.0
          0.0
          0.0]
    A = [10.0 -7.0 -1.0 0.0
          1.0  0.5  0.0 1.0]
    b = [5.0
         3.0]
    obj, X, y, Z, info, runhist = sdpt3(blk, [Matrix(A')], [c], b)
    tol = 1e-5
    @test obj[1] ≈ -1/8 atol=tol rtol=tol
    @test obj[2] ≈ -1/8 atol=tol rtol=tol
    @test X[1] ≈ [47/24, 25//12, 0, 0] atol=tol rtol=tol
    @test y ≈ [1/8, -1/4] atol=tol rtol=tol
    @test Z[1] ≈ [0, 0, 1/8, 1/4] atol=tol rtol=tol
    @test info["termcode"] == 0.0
end


text:
# This file was generated, do not modify it. # hide
train = 1:nrows(data_train)
test = last(train) .+ (1:nrows(data_test));


text:
using KerrQNMShifts
using Test

@testset "KerrQNMShifts.jl" begin
    # Write your tests here.
end

max_stars_repo_name: Potatoasad/KerrQNMShifts

max_issues_repo_name: Potatoasad/KerrQNMShifts

max_forks_repo_name: Potatoasad/KerrQNMShifts


text:
module Day19

using SparseArrays

_inp = [109,424,203,1,21102,11,1,0,1106,0,282,21101,18,0,0,1106,0,259,1202,1,1,221,203,1,21101,31,0,0,1105,1,282,21101,0,38,0,1106,0,259,20101,0,23,2,22101,0,1,3,21102,1,1,1,21102,57,1,0,1105,1,303,2102,1,1,222,21002,221,1,3,21001,221,0,2,21101,0,259,1,21101,80,0,0,1106,0,225,21101,0,149,2,21101,91,0,0,1106,0,303,1202,1,1,223,21002,222,1,4,21102,1,259,3,21101,225,0,2,21101,0,225,1,21101,118,0,0,1106,0,225,21001,222,0,3,21102,1,58,2,21102,1,133,0,1106,0,303,21202,1,-1,1,22001,223,1,1,21101,148,0,0,1105,1,259,1201,1,0,223,21002,221,1,4,20102,1,222,3,21101,21,0,2,1001,132,-2,224,1002,224,2,224,1001,224,3,224,1002,132,-1,132,1,224,132,224,21001,224,1,1,21102,195,1,0,105,1,109,20207,1,223,2,20102,1,23,1,21101,0,-1,3,21102,1,214,0,1106,0,303,22101,1,1,1,204,1,99,0,0,0,0,109,5,2101,0,-4,249,22101,0,-3,1,21201,-2,0,2,22101,0,-1,3,21101,250,0,0,1106,0,225,22101,0,1,-4,109,-5,2106,0,0,109,3,22107,0,-2,-1,21202,-1,2,-1,21201,-1,-1,-1,22202,-1,-2,-2,109,-3,2105,1,0,109,3,21207,-2,0,-1,1206,-1,294,104,0,99,22101,0,-2,-2,109,-3,2106,0,0,109,5,22207,-3,-4,-1,1206,-1,346,22201,-4,-3,-4,21202,-3,-1,-1,22201,-4,-1,2,21202,2,-1,-1,22201,-4,-1,1,22101,0,-2,3,21101,0,343,0,1105,1,303,1105,1,415,22207,-2,-3,-1,1206,-1,387,22201,-3,-2,-3,21202,-2,-1,-1,22201,-3,-1,3,21202,3,-1,-1,22201,-3,-1,2,22102,1,-4,1,21102,1,384,0,1105,1,303,1105,1,415,21202,-4,-1,-4,22201,-4,-3,-4,22202,-3,-2,-2,22202,-2,-4,-4,22202,-3,-2,-3,21202,-4,-1,-2,22201,-3,-2,1,21202,1,1,-4,109,-5,2105,1,0]

inp = Dict(zip(0:length(_inp) - 1, Ref.(_inp))) # convert to 0-based index

newref() = Ref(0)

function simulate!(v, inc, outc)
    i = 0
    relbase = 0

    function getval(modes, idx)
        modes[idx] == 0 && return get!(newref, v, v[i + idx][])
        modes[idx] == 1 && return get(v, i + idx, 0)
        modes[idx] == 2 && return get!(newref, v, relbase + v[i + idx][])
    end

    while true
        opcode, modes = v[i][] % 100, digits(div(v[i][], 100), pad = 3)

        function rtype!(op)
            getval(modes, 3)[] = op(getval(modes, 1)[], getval(modes, 2)[])
            i += 4
        end

        jtype!(op) = i = op(getval(modes, 1)[], 0) ? getval(modes, 2)[] : i + 3

        function itype!(op)
            getval(modes, 3)[] = Int(op(getval(modes, 1)[], getval(modes, 2)[]))
            i += 4
        end

        opcode in 1:9 || return
        opcode == 1 && rtype!(+)
        opcode == 2 && rtype!(*)
        opcode == 5 && jtype!(!=)
        opcode == 6 && jtype!(==)
        opcode == 7 && itype!(<)
        opcode == 8 && itype!(==)

        if opcode == 3
            getval(modes, 1)[] = take!(inc)
            i += 2
        elseif opcode == 4
            put!(outc, getval(modes, 1)[])
            i += 2
        elseif opcode == 9
            relbase += getval(modes, 1)[]
            i += 2
        end
    end
end

function solvea(v)
    inc, outc = Channel(2), Channel(1)

    cnt = 0
    for x in 0:49, y in 0:49
        put!(inc, x)
        put!(inc, y)
        simulate!(deepcopy(v), inc, outc)
        cnt += take!(outc)
    end

    cnt
end

function solveb(v) #incomplete
    inc, outc = Channel(2), Channel(1)

    x, y = 99900, 99900
    mat = spzeros(10000, 10000)
    seen = Set()
    reasons = []

    while true
        length(reasons) == 2 && return x * 10000 + y

        for i in x:x + 100, j in y:y + 100
            (i, j) in seen && continue
            t = Threads.@spawn simulate!(deepcopy(v), inc, outc)
            put!(inc, x)
            put!(inc, y)
            push!(seen, (i, j))
            mat[i, j] = take!(outc)

        end

    end

    for x in 0:601, y in 0:601
    end

    for x in 0:500, y in 0:500
        sum(mat[x:x + 100, y:y + 100]) == 10000 && return x * 10000 + y
    end
end

end

max_stars_repo_name: AnAverageHuman/competitive

max_issues_repo_name: AnAverageHuman/competitive

max_forks_repo_name: AnAverageHuman/competitive