ideas-statistics (empty) → 1.0
raw patch · 21 files changed
+4068/−0 lines, 21 filesdep +basedep +containersdep +ideassetup-changed
Dependencies added: base, containers, ideas, ideas-math-types
Files
- LICENSE.txt +202/−0
- NOTICE.txt +48/−0
- Setup.lhs +4/−0
- ideas-statistics.cabal +62/−0
- src/Domain/Hypothesis/BuggyRules.hs +406/−0
- src/Domain/Hypothesis/Common.hs +267/−0
- src/Domain/Hypothesis/Constraints.hs +434/−0
- src/Domain/Hypothesis/Examples.hs +450/−0
- src/Domain/Hypothesis/Exercises.hs +178/−0
- src/Domain/Hypothesis/Rules.hs +666/−0
- src/Domain/Hypothesis/Strategies.hs +74/−0
- src/Domain/Hypothesis/Tables.hs +88/−0
- src/Domain/Statistics/Component.hs +328/−0
- src/Domain/Statistics/ComponentSet.hs +163/−0
- src/Domain/Statistics/Data.hs +86/−0
- src/Domain/Statistics/Parser.hs +71/−0
- src/Domain/Statistics/Rules.hs +94/−0
- src/Domain/Statistics/Symbols.hs +133/−0
- src/Domain/Statistics/Views.hs +183/−0
- src/Main.hs +38/−0
- src/Service/ConstraintServices.hs +93/−0
+ LICENSE.txt view
@@ -0,0 +1,202 @@+ + Apache License + Version 2.0, January 2004 + http://www.apache.org/licenses/ + + TERMS AND CONDITIONS FOR USE, REPRODUCTION, AND DISTRIBUTION + + 1. Definitions. + + "License" shall mean the terms and conditions for use, reproduction, + and distribution as defined by Sections 1 through 9 of this document. + + "Licensor" shall mean the copyright owner or entity authorized by + the copyright owner that is granting the License. + + "Legal Entity" shall mean the union of the acting entity and all + other entities that control, are controlled by, or are under common + control with that entity. For the purposes of this definition, + "control" means (i) the power, direct or indirect, to cause the + direction or management of such entity, whether by contract or + otherwise, or (ii) ownership of fifty percent (50%) or more of the + outstanding shares, or (iii) beneficial ownership of such entity. + + "You" (or "Your") shall mean an individual or Legal Entity + exercising permissions granted by this License. + + "Source" form shall mean the preferred form for making modifications, + including but not limited to software source code, documentation + source, and configuration files. + + "Object" form shall mean any form resulting from mechanical + transformation or translation of a Source form, including but + not limited to compiled object code, generated documentation, + and conversions to other media types. + + "Work" shall mean the work of authorship, whether in Source or + Object form, made available under the License, as indicated by a + copyright notice that is included in or attached to the work + (an example is provided in the Appendix below). + + "Derivative Works" shall mean any work, whether in Source or Object + form, that is based on (or derived from) the Work and for which the + editorial revisions, annotations, elaborations, or other modifications + represent, as a whole, an original work of authorship. For the purposes + of this License, Derivative Works shall not include works that remain + separable from, or merely link (or bind by name) to the interfaces of, + the Work and Derivative Works thereof. + + "Contribution" shall mean any work of authorship, including + the original version of the Work and any modifications or additions + to that Work or Derivative Works thereof, that is intentionally + submitted to Licensor for inclusion in the Work by the copyright owner + or by an individual or Legal Entity authorized to submit on behalf of + the copyright owner. For the purposes of this definition, "submitted" + means any form of electronic, verbal, or written communication sent + to the Licensor or its representatives, including but not limited to + communication on electronic mailing lists, source code control systems, + and issue tracking systems that are managed by, or on behalf of, the + Licensor for the purpose of discussing and improving the Work, but + excluding communication that is conspicuously marked or otherwise + designated in writing by the copyright owner as "Not a Contribution." + + "Contributor" shall mean Licensor and any individual or Legal Entity + on behalf of whom a Contribution has been received by Licensor and + subsequently incorporated within the Work. + + 2. Grant of Copyright License. Subject to the terms and conditions of + this License, each Contributor hereby grants to You a perpetual, + worldwide, non-exclusive, no-charge, royalty-free, irrevocable + copyright license to reproduce, prepare Derivative Works of, + publicly display, publicly perform, sublicense, and distribute the + Work and such Derivative Works in Source or Object form. + + 3. Grant of Patent License. Subject to the terms and conditions of + this License, each Contributor hereby grants to You a perpetual, + worldwide, non-exclusive, no-charge, royalty-free, irrevocable + (except as stated in this section) patent license to make, have made, + use, offer to sell, sell, import, and otherwise transfer the Work, + where such license applies only to those patent claims licensable + by such Contributor that are necessarily infringed by their + Contribution(s) alone or by combination of their Contribution(s) + with the Work to which such Contribution(s) was submitted. If You + institute patent litigation against any entity (including a + cross-claim or counterclaim in a lawsuit) alleging that the Work + or a Contribution incorporated within the Work constitutes direct + or contributory patent infringement, then any patent licenses + granted to You under this License for that Work shall terminate + as of the date such litigation is filed. + + 4. Redistribution. You may reproduce and distribute copies of the + Work or Derivative Works thereof in any medium, with or without + modifications, and in Source or Object form, provided that You + meet the following conditions: + + (a) You must give any other recipients of the Work or + Derivative Works a copy of this License; and + + (b) You must cause any modified files to carry prominent notices + stating that You changed the files; and + + (c) You must retain, in the Source form of any Derivative Works + that You distribute, all copyright, patent, trademark, and + attribution notices from the Source form of the Work, + excluding those notices that do not pertain to any part of + the Derivative Works; and + + (d) If the Work includes a "NOTICE" text file as part of its + distribution, then any Derivative Works that You distribute must + include a readable copy of the attribution notices contained + within such NOTICE file, excluding those notices that do not + pertain to any part of the Derivative Works, in at least one + of the following places: within a NOTICE text file distributed + as part of the Derivative Works; within the Source form or + documentation, if provided along with the Derivative Works; or, + within a display generated by the Derivative Works, if and + wherever such third-party notices normally appear. The contents + of the NOTICE file are for informational purposes only and + do not modify the License. You may add Your own attribution + notices within Derivative Works that You distribute, alongside + or as an addendum to the NOTICE text from the Work, provided + that such additional attribution notices cannot be construed + as modifying the License. + + You may add Your own copyright statement to Your modifications and + may provide additional or different license terms and conditions + for use, reproduction, or distribution of Your modifications, or + for any such Derivative Works as a whole, provided Your use, + reproduction, and distribution of the Work otherwise complies with + the conditions stated in this License. + + 5. Submission of Contributions. Unless You explicitly state otherwise, + any Contribution intentionally submitted for inclusion in the Work + by You to the Licensor shall be under the terms and conditions of + this License, without any additional terms or conditions. + Notwithstanding the above, nothing herein shall supersede or modify + the terms of any separate license agreement you may have executed + with Licensor regarding such Contributions. + + 6. Trademarks. This License does not grant permission to use the trade + names, trademarks, service marks, or product names of the Licensor, + except as required for reasonable and customary use in describing the + origin of the Work and reproducing the content of the NOTICE file. + + 7. Disclaimer of Warranty. Unless required by applicable law or + agreed to in writing, Licensor provides the Work (and each + Contributor provides its Contributions) on an "AS IS" BASIS, + WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or + implied, including, without limitation, any warranties or conditions + of TITLE, NON-INFRINGEMENT, MERCHANTABILITY, or FITNESS FOR A + PARTICULAR PURPOSE. You are solely responsible for determining the + appropriateness of using or redistributing the Work and assume any + risks associated with Your exercise of permissions under this License. + + 8. Limitation of Liability. In no event and under no legal theory, + whether in tort (including negligence), contract, or otherwise, + unless required by applicable law (such as deliberate and grossly + negligent acts) or agreed to in writing, shall any Contributor be + liable to You for damages, including any direct, indirect, special, + incidental, or consequential damages of any character arising as a + result of this License or out of the use or inability to use the + Work (including but not limited to damages for loss of goodwill, + work stoppage, computer failure or malfunction, or any and all + other commercial damages or losses), even if such Contributor + has been advised of the possibility of such damages. + + 9. Accepting Warranty or Additional Liability. While redistributing + the Work or Derivative Works thereof, You may choose to offer, + and charge a fee for, acceptance of support, warranty, indemnity, + or other liability obligations and/or rights consistent with this + License. However, in accepting such obligations, You may act only + on Your own behalf and on Your sole responsibility, not on behalf + of any other Contributor, and only if You agree to indemnify, + defend, and hold each Contributor harmless for any liability + incurred by, or claims asserted against, such Contributor by reason + of your accepting any such warranty or additional liability. + + END OF TERMS AND CONDITIONS + + APPENDIX: How to apply the Apache License to your work. + + To apply the Apache License to your work, attach the following + boilerplate notice, with the fields enclosed by brackets "[]" + replaced with your own identifying information. (Don't include + the brackets!) The text should be enclosed in the appropriate + comment syntax for the file format. We also recommend that a + file or class name and description of purpose be included on the + same "printed page" as the copyright notice for easier + identification within third-party archives. + + Copyright [yyyy] [name of copyright owner] + + Licensed under the Apache License, Version 2.0 (the "License"); + you may not use this file except in compliance with the License. + You may obtain a copy of the License at + + http://www.apache.org/licenses/LICENSE-2.0 + + Unless required by applicable law or agreed to in writing, software + distributed under the License is distributed on an "AS IS" BASIS, + WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. + See the License for the specific language governing permissions and + limitations under the License.
+ NOTICE.txt view
@@ -0,0 +1,48 @@+Copyright 2020 Ideas project team++ Licensed under the Apache License, Version 2.0 (the "License");+ you may not use this file except in compliance with the License.+ You may obtain a copy of the License at++ http://www.apache.org/licenses/LICENSE-2.0++ Unless required by applicable law or agreed to in writing, software+ distributed under the License is distributed on an "AS IS" BASIS,+ WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.+ See the License for the specific language governing permissions and+ limitations under the License.+ +PACKAGE++ Interactive domain reasoner for statistics, based on the feedback services + of the Ideas framework. Provides feedback and hints on step-by-step + construction of hypothesis tests by students. Test types supported include + z-tests, t-tests for one group, dependent groups and independent groups, + correlation, ANOVA and Chi-square tests. The domain reasoner is used by the + learning environment Numworx (www.numworx.nl) and is available for + implementation in other learning environments as well. ++PROJECT++ Ideas (Interactive Domain-specific Exercise Assistants) is a joint research+ project between the Open University of the Netherlands and Utrecht+ University. The project's goal is to use software and compiler technology to+ build state-of-the-art components of intelligent tutoring systems (ITS) and+ learning environments. The 'ideas' software package provides a generic+ framework for constructing the expert knowledge module (also known as a+ domain reasoner) for an ITS or learning environment. More information about+ the project can be found on the project's homepage:+ + http://ideas.cs.uu.nl/+ +AUTHORS++ The authors and copyright holders of the 'ideas-statistics' software package are:+ - Bastiaan Heeren (bastiaan.heeren@ou.nl, maintainer)+ - Sietske Tacoma+ - Noeri Huisman++CREDITS++ Paul Drijvers, Johan Jeuring, Martijn Fleuren, Wim van Velthoven, Peter Boon,+ Jeltje Wassenberg-Severijnen, Corine Geurts, Alex Gerdes
+ Setup.lhs view
@@ -0,0 +1,4 @@+#! /usr/bin/env runhaskell++> import Distribution.Simple+> main = defaultMain
+ ideas-statistics.cabal view
@@ -0,0 +1,62 @@+name: ideas-statistics +version: 1.0 +synopsis: Interactive domain reasoner for statistics +homepage: http://ideas.cs.uu.nl/ +description: + + Interactive domain reasoner for statistics, based on the feedback services + of the Ideas framework. Provides feedback and hints on step-by-step construction + of hypothesis tests by students. Test types supported include z-tests, t-tests + for one group, dependent groups and independent groups, correlation, ANOVA and + Chi-square tests. The domain reasoner is used by the learning environment Numworx + (www.numworx.nl) and is available for implementation in other learning + environments as well. + +category: Education +copyright: (c) 2020 +license: Apache-2.0 +license-file: LICENSE.txt +author: Bastiaan Heeren, Sietske Tacoma, Noeri Huisman +maintainer: bastiaan.heeren@ou.nl +stability: provisional +extra-source-files: NOTICE.txt +build-type: Simple +cabal-version: >= 1.8.0.2 +tested-with: GHC == 8.6.5 + +source-repository head + type: svn + location: https://ideastest.science.uu.nl/svn/ideas/Tutors/statistics + +-------------------------------------------------------------------------------- + +executable statistics.cgi + main-is: Main.hs + hs-source-dirs: src + ghc-options: -Wall + + build-depends: + base >= 4.8 && < 5, + containers, + ideas >= 1.8, + ideas-math-types >= 1.1 + + other-modules: + Domain.Hypothesis.BuggyRules, + Domain.Hypothesis.Common, + Domain.Hypothesis.Constraints, + Domain.Hypothesis.Examples, + Domain.Hypothesis.Exercises, + Domain.Hypothesis.Rules, + Domain.Hypothesis.Strategies, + Domain.Hypothesis.Tables, + Domain.Statistics.Component, + Domain.Statistics.ComponentSet, + Domain.Statistics.Data, + Domain.Statistics.Parser, + Domain.Statistics.Rules, + Domain.Statistics.Symbols, + Domain.Statistics.Views, + Service.ConstraintServices + +--------------------------------------------------------------------------------
+ src/Domain/Hypothesis/BuggyRules.hs view
@@ -0,0 +1,406 @@+-----------------------------------------------------------------------------+-- Copyright 2020, Ideas project team. This file is distributed under the+-- terms of the Apache License 2.0. For more information, see the files+-- "LICENSE.txt" and "NOTICE.txt", which are included in the distribution.+-----------------------------------------------------------------------------+module Domain.Hypothesis.BuggyRules + ( buggyRules+ , buggyTValueSided, buggyRValueSided, buggyZValueSided+ , buggyTValueTestValue, buggyRValueTestValue, buggyZValueTestValue+ , buggyChiValueSided, buggyChiValueTestValue+ ) where++import Control.Monad+import Data.List+import Domain.Hypothesis.Common+import Domain.Math.Data.Relation+import Domain.Math.Expr hiding ((.*.), (./.))+import Domain.Math.Numeric.Views+import Domain.Hypothesis.Rules+import Domain.Statistics.ComponentSet+import Domain.Statistics.Views+import Ideas.Common.Library++----------------------------------------------------------+-- Buggy rules++-- Transform the steps described below into a rule, the pattern is factored out.+-- The list can be freely appended with new buggy rules+buggyRules:: [Rule ComponentSet]+buggyRules =+ let _pat (f, name, descr) = (describe descr . buggyRule name) f+ _name_common x = "buggy.component." ++ x+ in map _pat+ [ ( wrongHARule+ , _name_common "ha"+ , "Choosing wrong alternative hypothesis"+ )+ , ( wrongAlphaRule+ , _name_common "alpha"+ , "Buggy rule for adding the wrong alpha"+ )+ , ( wrongDf+ , _name_common "df"+ , "Buggy rule for adding the wrong degrees of freedom"+ )+ ] ++ [buggyHypothesesSampleMean, buggyHASampleMean, buggyH0SampleMean]+ ++ buggyRulesTest+ ++ [ buggyTestZValueRule, buggyTestTValueRule+ , buggyTValuePositive+ , buggyZValueAlpha, buggyZValueSided, buggyZValueTestValue+ , buggyTValueAlpha, buggyTValueSided, buggyTValueDf, buggyTValueTestValue+ , buggyRValueAlpha, buggyRValueSided, buggyRValueDf, buggyRValueTestValue+ , buggyFValueAlpha, buggyFValueDf, buggyFValueTestValue+ , buggyChiValueAlpha, buggyChiValueSided, buggyChiValueDf, buggyChiValueTestValue+ , buggyTestValueCritical+ ]++alphaValues :: [Double]+alphaValues = [0.01, 0.05, 0.10]++alphaChiValues :: [Double]+alphaChiValues = [0.10, 0.05, 0.025, 0.01, 0.005]++sidedValues :: [Sided]+sidedValues = [TwoSided, LeftSided, RightSided]++incorrectDf :: Double -> [Double]+incorrectDf correctDf = [correctDf - 1, correctDf + 1, correctDf + 2]++confusingTestValue :: ComponentSet -> [Double]+confusingTestValue cs = do + r <- inferTestValue cs + matchM doubleView (rightHandSide r)++buggyTestValueCritical :: Rule ComponentSet+buggyTestValueCritical =+ describe "Buggy rule for confusing the test value with critical value" $+ buggyRule "buggy.component.test-value.critical" f+ where+ f :: ComponentSet -> Maybe ComponentSet+ f cs = do+ guard (derived cs `doesNotContain` TestValue)+ test <- inferTestChoice cs+ var <- inferVar cs+ alpha <- getExpr SignificanceLevel cs >>= match doubleView+ crit <- + case test of+ Anova -> do+ (dfBetween, dfWithin) <- inferDfBetweenWithin cs+ computeCriticalF dfBetween dfWithin alpha+ ZTest -> do+ sided <- inferSidedness cs+ return $ computeCriticalZ sided alpha+ RPearson -> do+ sided <- inferSidedness cs+ df <- inferDf cs >>= match doubleView+ return $ computeCriticalR sided alpha df+ _ | isTTest test -> do+ sided <- inferSidedness cs+ df <- inferDf cs >>= match doubleView+ return $ computeCriticalT sided alpha df+ _ -> Nothing+ return $ append TestValue (CRelation $ var .==. toExpr crit) cs++------------------------------------------------------++buggyZValueAlpha :: Rule ComponentSet+buggyZValueAlpha =+ describe "Buggy rule for looking up a z-value with wrong alpha" .+ buggyRule "buggy.component.critical.z-value.alpha" $+ inferCriticalZWith $ \sided alpha ->+ [ computeCriticalZ sided a | a <- alphaValues, a /= alpha ]++buggyZValueSided :: Rule ComponentSet+buggyZValueSided =+ describe "Buggy rule for looking up a z-value with wrong sidedness" .+ buggyRule "buggy.component.critical.z-value.sided" $+ inferCriticalZWith $ \sided alpha ->+ [ computeCriticalZ s alpha | s <- sidedValues, s /= sided ]++buggyZValueTestValue :: Rule ComponentSet+buggyZValueTestValue =+ describe "Buggy rule for confusing z-value with test-value" .+ buggyRule "buggy.component.critical.z-value.test-value" $ \cs -> + inferCriticalZWith (\_ _ -> confusingTestValue cs) cs++------------------------------------------------------++buggyTValueAlpha :: Rule ComponentSet+buggyTValueAlpha =+ describe "Buggy rule for looking up a t-value with wrong alpha" .+ buggyRule "buggy.component.critical.t-value.alpha" $+ inferCriticalTWith $ \sided alpha df ->+ [ computeCriticalT sided a df | a <- alphaValues, a /= alpha ]++buggyTValueSided :: Rule ComponentSet+buggyTValueSided =+ describe "Buggy rule for looking up a t-value with wrong sidedness" .+ buggyRule "buggy.component.critical.t-value.sided" $+ inferCriticalTWith $ \sided alpha df ->+ [ computeCriticalT s alpha df | s <- sidedValues, s /= sided ]++buggyTValueDf :: Rule ComponentSet+buggyTValueDf =+ describe "Buggy rule for looking up a t-value with wrong df" .+ buggyRule "buggy.component.critical.t-value.df" $+ inferCriticalTWith $ \sided alpha df ->+ [ computeCriticalT sided alpha dfWrong | dfWrong <- incorrectDf df ]++buggyTValueTestValue :: Rule ComponentSet+buggyTValueTestValue =+ describe "Buggy rule for confusing t-value with test-value" .+ buggyRule "buggy.component.critical.t-value.test-value" $ \cs -> + inferCriticalTWith (\_ _ _ -> confusingTestValue cs) cs++------------------------------------------------------++buggyRValueAlpha :: Rule ComponentSet+buggyRValueAlpha =+ describe "Buggy rule for looking up a r-value with wrong alpha" .+ buggyRule "buggy.component.critical.r-value.alpha" $+ inferCriticalRWith $ \sided alpha df ->+ [ computeCriticalR sided a df | a <- alphaValues, a /= alpha ]++buggyRValueSided :: Rule ComponentSet+buggyRValueSided =+ describe "Buggy rule for looking up a r-value with wrong sidedness" .+ buggyRule "buggy.component.critical.r-value.sided" $+ inferCriticalRWith $ \sided alpha df ->+ [ computeCriticalR s alpha df | s <- sidedValues, s /= sided ]++buggyRValueDf :: Rule ComponentSet+buggyRValueDf =+ describe "Buggy rule for looking up a r-value with wrong df" .+ buggyRule "buggy.component.critical.r-value.df" $+ inferCriticalRWith $ \sided alpha df ->+ [ computeCriticalR sided alpha dfWrong | dfWrong <- incorrectDf df ]++buggyRValueTestValue :: Rule ComponentSet+buggyRValueTestValue =+ describe "Buggy rule for confusing r-value with test-value" .+ buggyRule "buggy.component.critical.r-value.test-value" $ \cs -> + inferCriticalRWith (\_ _ _ -> confusingTestValue cs) cs++------------------------------------------------------++buggyFValueAlpha :: Rule ComponentSet+buggyFValueAlpha =+ describe "Buggy rule for looking up f-value with wrong alpha" .+ buggyRule "buggy.component.critical.f-value.alpha" $+ inferCriticalFWith $ \dfBetween dfWithin alpha ->+ [ result+ | a <- [0.01, 0.05], a /= alpha + , result <- computeCriticalF dfBetween dfWithin a + ]++buggyFValueDf :: Rule ComponentSet+buggyFValueDf =+ describe "Buggy rule for looking up f-value with wrong df" .+ buggyRule "buggy.component.critical.f-value.df" $+ inferCriticalFWith $ \dfBetween dfWithin alpha ->+ [ result+ | dfB <- [dfBetween, 2]+ , dfW <- [dfWithin, 40]+ , dfBetween /= dfB || dfWithin /= dfW+ , result <- computeCriticalF dfB dfW alpha+ ]++buggyFValueTestValue :: Rule ComponentSet+buggyFValueTestValue =+ describe "Buggy rule for confusing f-value with test-value" .+ buggyRule "buggy.component.critical.f-value.test-value" $ \cs -> + inferCriticalFWith (\_ _ _ -> confusingTestValue cs) cs++------------------------------------------------------++buggyChiValueSided :: Rule ComponentSet+buggyChiValueSided =+ describe "Buggy rule for looking up a chi-value with wrong sidedness" .+ buggyRule "buggy.component.critical.chi-value.sided" $+ inferCriticalChiWith $ \sided alpha df ->+ [ result+ | s <- sidedValues+ , s /= sided + , result <- computeCriticalChi s alpha df + ]++buggyChiValueAlpha :: Rule ComponentSet+buggyChiValueAlpha =+ describe "Buggy rule for looking up chi-value with wrong alpha" .+ buggyRule "buggy.component.critical.chi-value.alpha" $+ inferCriticalChiWith $ \sided alpha df ->+ [ result+ | a <- alphaChiValues, a /= alpha + , result <- computeCriticalChi sided a df+ ]++buggyChiValueDf :: Rule ComponentSet+buggyChiValueDf =+ describe "Buggy rule for looking up chi-value with wrong df" .+ buggyRule "buggy.component.critical.chi-value.df" $+ inferCriticalChiWith $ \sided alpha df ->+ [ result + | dfWrong <- incorrectDf df + , result <- computeCriticalChi sided alpha dfWrong + ]++buggyChiValueTestValue :: Rule ComponentSet+buggyChiValueTestValue =+ describe "Buggy rule for confusing chi-value with test-value" .+ buggyRule "buggy.component.critical.chi-value.test-value" $ \cs -> + inferCriticalChiWith (\_ _ _ -> confusingTestValue cs) cs++------------------------------------------------------++buggyHypothesesSampleMean :: Rule ComponentSet+buggyHypothesesSampleMean =+ describe "Buggy rule for adding one or two hypotheses using sample mean instead of population mean" $ + buggyRule "buggy.component.hypotheses-samplemean" f+ where+ f :: ComponentSet -> [ComponentSet]+ f cs = do+ guard (derived cs `doesNotContain` NullHypothesis)+ guard (derived cs `doesNotContain` AlternativeHypothesis)+ ha <- getRelation AlternativeHypothesis cs+ let h0 = h0FromHA ha+ sm <- getExpr SampleMean cs+ rtp <- nub [EqualTo, relationType h0]+ let h0' = makeType rtp (leftHandSide h0) (rightHandSide h0)+ buggyH0 = makeType rtp (leftHandSide h0) sm+ buggyHA = makeType (relationType ha) (leftHandSide ha) sm+ results = + [ append AlternativeHypothesis (CRelation buggyHA) $+ append NullHypothesis (CRelation buggyH0) cs+ , append AlternativeHypothesis (CRelation ha) $+ append NullHypothesis (CRelation buggyH0) cs+ , append AlternativeHypothesis (CRelation buggyHA) $+ append NullHypothesis (CRelation h0') cs+ ]+ results++buggyH0SampleMean :: Rule ComponentSet+buggyH0SampleMean =+ describe "Buggy rule for adding h0 using sample mean instead of population mean" $ + buggyRule "buggy.component.h0-samplemean" f+ where+ f :: ComponentSet -> [ComponentSet]+ f cs = do+ guard (derived cs `doesNotContain` NullHypothesis)+ ha <- getRelation AlternativeHypothesis cs+ let h0 = h0FromHA ha+ sm <- getExpr SampleMean cs+ rtp <- nub [EqualTo, relationType h0]+ let buggyH0 = makeType rtp (leftHandSide h0) sm+ return $ append NullHypothesis (CRelation buggyH0) cs++buggyHASampleMean :: Rule ComponentSet+buggyHASampleMean =+ describe "Buggy rule for adding ha using sample mean instead of population mean" $ + buggyRule "buggy.component.ha-samplemean" f+ where+ f :: ComponentSet -> Maybe ComponentSet+ f cs = do+ guard (derived cs `doesNotContain` AlternativeHypothesis)+ ha <- getRelation AlternativeHypothesis cs+ sm <- getExpr SampleMean cs+ let buggyHA = makeType (relationType ha) (leftHandSide ha) sm+ return $ append AlternativeHypothesis (CRelation buggyHA) cs++buggyTestZValueRule :: Rule ComponentSet+buggyTestZValueRule =+ describe "Standard error instead of standard deviation" $+ buggyRule "buggy.component.test-z-value" f+ where+ f :: ComponentSet -> Maybe ComponentSet+ f cs = do+ guard (cs `doesNotContain` TestValue)+ let expr = (Var "M" - Var "mu") / Var "sigma"+ val <- match doubleView (getSubstitution cs |-> expr)+ var <- inferVar cs+ return $ append TestValue (CRelation (var .==. toExpr val)) cs++-- see buggyTestZValueRule: only difference is 's' instead of 'sigma'+buggyTestTValueRule :: Rule ComponentSet+buggyTestTValueRule =+ describe "Standard error instead of standard deviation" $ + buggyRule "buggy.component.test-t-value" f+ where+ f :: ComponentSet -> Maybe ComponentSet+ f cs = do+ guard (cs `doesNotContain` TestValue)+ let expr = (Var "M" - Var "mu") / Var "s"+ val <- match doubleView (getSubstitution cs |-> expr)+ var <- inferVar cs+ return $ append TestValue (CRelation (var .==. toExpr val)) cs++buggyTValuePositive :: Rule ComponentSet+buggyTValuePositive =+ buggyRule "buggy.component.critical.t-value-positive" f+ where+ f :: ComponentSet -> Maybe ComponentSet+ f cs = do+ cs' <- apply lookupTValueRule cs+ rel <- getRelation Critical cs'+ val <- match doubleView (rightHandSide rel)+ guard (val < 0)+ return $ append Critical (CRelation $ Var "tcrit" .==. toExpr (abs val)) cs++buggyRulesTest :: [Rule ComponentSet]+buggyRulesTest =+ [ buggyRule (makeId wrong correct) (wrongTest wrong correct)+ | wrong <- types+ , correct <- types+ , wrong /= correct+ ]+ where+ types = [TTestOne, TTestTwo, TTestPaired, ZTest]+ makeId t1 t2 = "buggy.component." ++ show t1 ++ ".should-be-" ++ show t2++wrongTest :: TestType -> TestType -> ComponentSet -> Maybe ComponentSet+wrongTest wrong correct cs = do+ guard (derived cs `doesNotContain` TestChoice)+ let tests = validTests cs+ guard (correct `elem` tests)+ return $ append TestChoice (CChoice $ TestType wrong) cs++-- Wrong sidedness in alternative hypothesis (due to sample mean)+wrongHARule :: ComponentSet -> Maybe ComponentSet+wrongHARule cs =+ let relFromMean :: Double -> Double -> RelationType+ relFromMean mean mu+ | mean < mu = LessThanOrEqualTo+ | otherwise = GreaterThanOrEqualTo+ in do+ guard (cs `doesNotContain` AlternativeHypothesis)+ ha <- getRelation AlternativeHypothesis cs -- used to be Claim+ let h0 = h0FromHA ha+ h0Rel = relationType h0+ guard (h0Rel == EqualTo)+ let lhs = leftHandSide ha+ rhs = rightHandSide ha+ rhs' <- match doubleView rhs+ mean <- match doubleView =<< getExpr SampleMean cs+ let rel = relFromMean mean rhs'+ return . append AlternativeHypothesis+ (CRelation (makeType rel lhs rhs)) $ cs++wrongAlphaRule :: ComponentSet -> Maybe ComponentSet+wrongAlphaRule cs = do+ guard (cs `doesNotContain` SignificanceLevel)+ return $ append SignificanceLevel (CExpr $ Number 0.1) cs++-- | Use the wrong df: t-tests are swapped so n instead of n-1 or n_1 + n_2 - 1+-- instead of n_1 + n_2 - 2+wrongDf :: ComponentSet -> Maybe ComponentSet+wrongDf cs =+ let dfWrong :: TestType -> Maybe Expr+ dfWrong TTestOne = Just (Var "n")+ dfWrong TTestPaired = Just (Var "n")+ dfWrong TTestTwo = Just (Var "n1" + Var "n2" - 1)+ dfWrong _ = Nothing+ in do+ guard (cs `doesNotContain` Df)+ test <- inferTestChoice cs+ df <- dfWrong test+ return $ append Df (CExpr df) cs
+ src/Domain/Hypothesis/Common.hs view
@@ -0,0 +1,267 @@+-----------------------------------------------------------------------------+-- Copyright 2020, Ideas project team. This file is distributed under the+-- terms of the Apache License 2.0. For more information, see the files+-- "LICENSE.txt" and "NOTICE.txt", which are included in the distribution.+-----------------------------------------------------------------------------+module Domain.Hypothesis.Common+ ( TypeOfTest(..),+ dependent, independent, chiSquared,+ pickAlpha, sidedFromHA,+ h0FromHA, h0FromHAEqualSign,+ validTests, testFormulaFromTest,+ degreesOfFreedomFromTest,+ computeCritical, computeCriticalT, computeCriticalZ+ , computeCriticalR, computeCriticalF, computeCriticalChi,+ computePValue, computePValueT, computePValueZ,+ chooseTypeOfTest, isTTest+ ) where++import Data.Maybe+import Ideas.Common.Library (Term(TCon))+import Domain.Statistics.Symbols+import Domain.Hypothesis.Tables+import Domain.Math.Expr hiding ((^))+import Domain.Math.Data.Relation+import Domain.Statistics.Component+import Domain.Statistics.ComponentSet+import Ideas.Common.View+import Domain.Math.Numeric.Views++data TypeOfTest = TestMean+ | CompareMeans+ | CompareMeansPaired+ | TestProportion+ | CompareProportions+ | TestCorrelation+ deriving (Eq, Show)++dependent, independent, chiSquared :: Expr+dependent = toExpr (TCon dependentSymbol [])+independent = toExpr (TCon independentSymbol [])+chiSquared = Var "chisq"++-- Pick a value for alpha+pickAlpha :: Double+pickAlpha = 0.05++-- Determine the sidedness based on the alternative hypothesis.+sidedFromHA :: Relation Expr -> Sided+sidedFromHA ha = fromRelation (relationType ha)+ where+ fromRelation :: RelationType -> Sided+ fromRelation LessThan = LeftSided+ fromRelation LessThanOrEqualTo = LeftSided+ fromRelation GreaterThan = RightSided+ fromRelation GreaterThanOrEqualTo = RightSided+ fromRelation EqualTo = TwoSided+ fromRelation NotEqualTo = TwoSided+ fromRelation rt = error $ "sidedFromHA: " ++ show rt++h0FromHA :: Relation Expr -> Relation Expr+h0FromHA ha = let h0Rel = inverseRelType (relationType ha)+ in makeType h0Rel (leftHandSide ha) (rightHandSide ha)++h0FromHAEqualSign :: Relation Expr -> Relation Expr+h0FromHAEqualSign ha = makeType EqualTo (leftHandSide ha) (rightHandSide ha)+ +inverseRelType :: RelationType -> RelationType+inverseRelType relType = fromMaybe relType (lookup relType table)+ where+ table = pairs ++ map (\(a,b) -> (b,a)) pairs+ pairs = [(LessThan, GreaterThanOrEqualTo), (LessThanOrEqualTo, GreaterThan), (EqualTo, NotEqualTo)]++-- Returns the valid tests for this state+-- Corner case for two samples (n1,n2 instead of n)+validTests :: ComponentSet -> [TestType]+validTests cs = + case getTestType TestChoice (initials cs) of+ Just x -> [x]+ _ -> select (chooseTypeOfTest cs)+ where+ n = fromMaybe 0 $ do + expr <- getExpr SampleSize cs+ match naturalView expr+ + psdKnown = cs `contains` PopulationSdev+ + select TestMean + | psdKnown = [ZTest]+ | n > 100 = [ZTest] -- Sietske: ignore TTestOne for now+ | otherwise = [TTestOne]+ select CompareMeans = [TTestTwo]+ select CompareMeansPaired = [TTestPaired]+ select TestProportion = [ZTest]+ select CompareProportions = [ZTest]+ select TestCorrelation = [RPearson, TTestOne]++-- Returns the test statistic formula from the chose test.+-- Note: this only works for the case of testing the mean, other cases should +-- give the formula as part of the exercise+testFormulaFromTest :: Monad m => TestType -> TypeOfTest -> m Expr+testFormulaFromTest testType typeOfTest = + case (testType, typeOfTest) of + (TTestOne, TestMean) -> + return $ (Var "M" - Var "mu") / (Var "s" / sqrt (Var "n"))+ (ZTest, TestMean) -> + return $ (Var "M" - Var "mu") / (Var "sigma" / sqrt (Var "n"))+ (TTestTwo, CompareMeans) -> + return $ (mean1 - mean2) / sqrt (toExpr PooledVariance * (1 / Var "n1" + 1 / Var "n2"))+ (TTestPaired, TestMean) ->+ return $ (Var "M" - Var "mu") / (Var "s" / sqrt (Var "n"))+ (TTestPaired, CompareMeansPaired) -> + return $ (mean1 - mean2) / (toExpr SampleSdev * sqrt (1 / Var "n")) -- to do: ask Sietske+ (_, TestProportion) -> + return $ (Var "p" - Var "p0") / sqrt (Var "p0" * (1.0 - Var "p0") / Var "n")+ (_, CompareProportions) -> + return $ (Var "p1" - Var "p2" - Var "d0") / sqrt (Var "p0" * (1.0 - Var "p0") / (Var "n" / 2))+ (TTestOne, TestCorrelation) -> + return $ (Var "r" * sqrt (Var "n" - 2)) / sqrt (1 - Var "r"**2)+ (RPearson, TestCorrelation) -> + return $ Var "r"+ _ -> + fail $ "teststatisticFromTest " ++ show (testType, typeOfTest)+ where+ mean1 = toExpr (One SampleMean)+ mean2 = toExpr (Two SampleMean)++-- Returns the formula for the degrees of degreesOfFreedom+degreesOfFreedomFromTest :: Monad m => TestType -> TypeOfTest -> m Expr+degreesOfFreedomFromTest TTestOne TestCorrelation = return $ Var "n" - 2+degreesOfFreedomFromTest TTestOne _ = return $ Var "n" - 1+degreesOfFreedomFromTest TTestPaired _ = return $ Var "n" - 1+degreesOfFreedomFromTest TTestTwo _ = return $ Var "n1" + Var "n2" - 2+degreesOfFreedomFromTest RPearson _ = return $ Var "n" - 2+degreesOfFreedomFromTest _ _ = fail "degrees of freedom test failed"++-- Returns the critical value from the given test and alpha+computeCritical :: Monad m => TestType -> Sided -> Double -> Maybe Double -> m Double+computeCritical test sided alpha mdf+ | isTTest test = case mdf of+ Just df -> return $ computeCriticalT sided alpha df+ Nothing -> fail "df missing"+ | test == RPearson = case mdf of+ Just df -> return $ computeCriticalR sided alpha df+ Nothing -> fail "df missing"+ | test == ZTest = return $ computeCriticalZ sided alpha+ | otherwise = fail "unknown test"++isTTest :: TestType -> Bool+isTTest TTestOne = True+isTTest TTestTwo = True+isTTest TTestPaired = True+isTTest _ = False++computeCriticalR :: Sided -> Double -> Double -> Double+computeCriticalR sided alpha df = sqrt (t ^ (2 :: Int) / (t ^ (2 :: Int) + df)) + where+ t = computeCriticalT sided alpha df++computeCriticalF :: Monad m => Double -> Double -> Double -> m Double+computeCriticalF dfBetween dfWithin alpha = + maybe (fail "unknown critical-f value") return (fTable dfBetween dfWithin alpha) -- TODO Sietske+++computeCriticalChi :: Monad m => Sided -> Double -> Double -> m Double+computeCriticalChi TwoSided alpha df = chivalue' df (alpha / 2)+computeCriticalChi LeftSided alpha df = negate <$> chivalue' df alpha+computeCriticalChi RightSided alpha df = chivalue' df alpha++chivalue' :: Monad m => Double -> Double -> m Double+chivalue' df alpha = + maybe (fail "unknown critical-chi value") return $ chiTable alpha (round df)++computeCriticalT :: Sided -> Double -> Double -> Double+computeCriticalT TwoSided alpha df = tvalue' df (alpha / 2)+computeCriticalT LeftSided alpha df = - tvalue' df alpha+computeCriticalT RightSided alpha df = tvalue' df alpha++tvalue' :: Double -> Double -> Double+tvalue' df alpha | isJust tableLookup = fromJust tableLookup+ | otherwise = fromInteger (round $ findValue (tvalue df) 0.00005 (0.5 - alpha) * 1000) / 1000.0+ where+ tableLookup = tTable alpha (round df)++computeCriticalZ :: Sided -> Double -> Double+computeCriticalZ TwoSided alpha = zvalue' (alpha / 2)+computeCriticalZ LeftSided alpha = - zvalue' alpha+computeCriticalZ RightSided alpha = zvalue' alpha++zvalue' :: Double -> Double+zvalue' alpha | isJust tableLookup = fromJust tableLookup+ | otherwise = fromInteger (round $ findValue zvalue 0.00005 (0.5 - alpha) * 1000) / 1000.0+ where+ tableLookup = zTable alpha++-- | Utils for computing t-/z-/p-values+zvalue :: Double -> Double+zvalue x = (1.0 / sqrt (2.0 * pi)) * exp 1 ** negate (x**2/2)++tvalue :: Double -> (Double -> Double)+tvalue df x = 1 / (sqrt df * beta 0.5 (df / 2.0)) * (1.0 + (x*x) / df) ** negate ((df + 1.0) / 2.0)++-- Utility functions for finding the t-value+-- Source: https://wiki.haskell.org/index.php?title=Gamma_and_Beta_function+cof :: [Double]+cof = [76.18009172947146,-86.50532032941677,24.01409824083091,-1.231739572450155,0.001208650973866179,-0.000005395239384953]+ +ser :: Double+ser = 1.000000000190015+ +gammaln :: Double -> Double+gammaln xx = let tmp' = (xx+5.5) - (xx+0.5)*log(xx+5.5)+ ser' = foldl (+) ser $ map (\(y,c) -> c/(xx+y)) $ zip [1..] cof+ in -tmp' + log(2.5066282746310005 * ser' / xx)+ +beta :: Double -> Double -> Double+beta z w = exp (gammaln z + gammaln w - gammaln (z+w))++findValue :: (Double -> Double) -> Double -> Double -> Double+findValue f stepSize target = fst $ until (\(_, x) -> x >= target) (\(a, x) -> (a + stepSize, x + stepSize * f a)) (0, 0)++findValue' :: (Double -> Double) -> Double -> Double -> Double+findValue' f stepSize target = snd $ until (\(a, _) -> a >= target) (\(a, x) -> (a + stepSize, x + stepSize * f a)) (0, 0)+++-- Compute the P-value+computePValue :: Monad m => TestType -> Sided -> Double -> Maybe Double -> m Double+computePValue test sided ts (Just df) + | isTTest test = return $ computePValueT sided ts df+ | otherwise = fail "cannot compute p-value"+computePValue ZTest sided ts _ = + return $ computePValueZ sided ts+computePValue _ _ _ _ = + fail "cannot compute p-value"++computePValueT :: Sided -> Double -> Double -> Double+computePValueT TwoSided ts df = (0.5 - pvalueT df (abs ts)) * 2+computePValueT LeftSided ts df | ts < 0 = 0.5 - pvalueT df (abs ts)+ | otherwise = 0.5 + pvalueT df ts+computePValueT RightSided ts df | ts < 0 = 1.0 - (0.5 - pvalueT df (abs ts))+ | otherwise = 0.5 - pvalueT df ts+ +pvalueT :: Double -> Double -> Double+pvalueT df testStatistic | testStatistic > 4.0 = 0.5 -- Performance gain, assuming precision is not needed for exceptionally high values+ | otherwise = findValue' (tvalue df) 0.00005 testStatistic++computePValueZ :: Sided -> Double -> Double+computePValueZ TwoSided ts = (0.5 - pvalueZ (abs ts)) * 2+computePValueZ LeftSided ts | ts < 0 = 0.5 - pvalueZ (abs ts)+ | otherwise = 0.5 + pvalueZ ts+computePValueZ RightSided ts | ts < 0 = 1.0 - (0.5 - pvalueZ (abs ts))+ | otherwise = 0.5 - pvalueZ ts++pvalueZ :: Double -> Double+pvalueZ testStatistic | testStatistic > 4.0 = 0.5 -- Performance gain, assuming precision is not needed for exceptionally high values+ | otherwise = findValue' zvalue 0.00005 testStatistic++-- | Utility function+chooseTypeOfTest :: ComponentSet -> TypeOfTest+chooseTypeOfTest cs + | contains cs Correlation = TestCorrelation+ | contains cs (Two SampleMean) = if contains cs (Two SampleSize)+ then CompareMeans+ else CompareMeansPaired+ | contains cs (Two Proportion) = CompareProportions+ | contains cs Proportion = TestProportion+ | contains cs SampleMean = TestMean+ | otherwise = TestMean
+ src/Domain/Hypothesis/Constraints.hs view
@@ -0,0 +1,434 @@+----------------------------------------------------------------------------- +-- Copyright 2020, Ideas project team. This file is distributed under the +-- terms of the Apache License 2.0. For more information, see the files +-- "LICENSE.txt" and "NOTICE.txt", which are included in the distribution. +----------------------------------------------------------------------------- +-- to do: split constraints in check against initial value, and inferred value +module Domain.Hypothesis.Constraints + ( hypothesisConstraints, equalDouble + ) where + +import Control.Monad +import Control.Applicative +import Data.Maybe +import Domain.Hypothesis.Common +import Domain.Hypothesis.Rules +import Domain.Math.Data.Relation +import Domain.Math.Numeric.Views +import Domain.Statistics.ComponentSet +import Ideas.Text.OpenMath.Dictionary.Arith1 +import Ideas.Common.Constraint +import Domain.Statistics.Views +import Ideas.Common.Library hiding (Predicate) +import Domain.Math.Expr.Data +import Prelude hiding (until) + +hypothesisConstraints :: [Constraint ComponentSet] +hypothesisConstraints = + [ checkHA + , checkH0 + , checkAlpha + , checkSided -- check this before checkCritical + , checkTest -- check this before checkCritical + , checkDf -- check this before checkCritical + , checkDfAnova -- check this before checkCritical + , checkGroups + , checkCritical, checkRejectionCritical + , checkSampleMean, checkSampleSdev, checkPopulationSdev + , checkPopulationMean, checkSampleSize, checkCorrelation + , checkStandardError + , checkTestFormula, checkTestValue + , checkPValue + -- checking conclusions should be last + , checkConclusionPValue, checkConclusionCritical, checkConclusionHypotheses + , checkObservedTotals, checkExpectedFrequencies + ] + +checkSampleMean, checkSampleSdev, checkPopulationSdev, + checkPopulationMean, checkSampleSize, checkCorrelation :: Constraint ComponentSet +checkSampleMean = checkForInitial "samplemean" SampleMean +checkSampleSdev = checkForInitial "samplesdev" SampleSdev +checkPopulationSdev = checkForInitial "populationsdev" PopulationSdev +checkPopulationMean = checkForInitial "populationmean" PopulationMean +checkSampleSize = checkForInitial "samplesize" SampleSize +checkCorrelation = checkForInitial "correlation" Correlation + +checkForInitial :: String -> ComponentId -> Constraint ComponentSet +checkForInitial s cid = makeConstraint ("check" # s) $ \cs -> do + relevance $ guard $ isInitial cid cs && isDerived cid cs + x <- get cid (initials cs) + y <- get cid (derived cs) + guard $ x == y + +checkAlpha :: Constraint ComponentSet +checkAlpha = makeConstraint "check.alpha" $ \cs -> do + relevance $ guard $ isDerived SignificanceLevel cs + if isInitial SignificanceLevel cs then do + -- initial alpha overrules the default value 0.05 + x <- get SignificanceLevel (initials cs) + y <- get SignificanceLevel (derived cs) + guard $ x == y + else do + -- use the default value 0.05 + alpha <- getExpr SignificanceLevel cs + val <- matchM doubleView alpha + guard $ val == 0.05 + +checkH0 :: Constraint ComponentSet +checkH0 = makeConstraint "check.h0" $ \cs -> do + relevance $ guard $ isDerived NullHypothesis cs + test <- inferTestChoice cs + case test of + ChiSquared -> do + h0 <- getExpr NullHypothesis cs + guard (h0 == independent) + _ -> do + h0 <- getRelation NullHypothesis cs + h0' <- h0FromHA <$> getRelation AlternativeHypothesis (initials cs) + unless (leftHandSide h0 == leftHandSide h0') $ + fail "parameter mismatch" + unless (rightHandSide h0 == rightHandSide h0') $ + fail "value mismatch" + unless (relationType h0 == relationType h0' || relationType h0 == EqualTo) $ + fail "sign mismatch" + +checkHA :: Constraint ComponentSet +checkHA = makeConstraint "check.ha" $ \cs -> do + relevance $ guard $ isDerived AlternativeHypothesis cs + test <- inferTestChoice cs + case test of + ChiSquared -> do + ha <- getExpr AlternativeHypothesis cs + guard (ha == dependent) + _ -> do + ha <- getRelation AlternativeHypothesis cs + ha' <- case get AlternativeHypothesis (initials cs) of + Just (CRelation x) -> return x + _ -> empty + unless (leftHandSide ha == leftHandSide ha') $ + fail "parameter mismatch" + unless (rightHandSide ha == rightHandSide ha') $ + fail "value mismatch" + unless (relationType ha == relationType ha') $ + fail "sign mismatch" + +checkConclusionPValue :: Constraint ComponentSet +checkConclusionPValue = makeConstraint "check.conclusion.p-value" $ \cs -> do + relevance $ guard $ isDerived ConclusionPValue cs + conc <- getRelation ConclusionPValue cs + let cs' = substitute cs + p <- matchM doubleView =<< getExpr PValue cs' + alpha <- matchM doubleView =<< getExpr SignificanceLevel cs' + guard $ if p <= alpha + then relationType conc == LessThanOrEqualTo + else relationType conc == GreaterThan + +checkSided :: Constraint ComponentSet +checkSided = makeConstraint "check.sided" $ \cs -> do + relevance $ guard $ isDerived Sidedness cs + -- TODO: make order of checking constraints explicit + -- && isNothing (isViolated checkHA cs) + sided <- getSided Sidedness cs + case getTestType TestChoice (initials cs) of + Just Anova -> + unless (sided == RightSided) (fail "anova") + Just ChiSquared -> + unless (sided == RightSided) (fail "chi-squared") + _ -> do + ha <- getRelation AlternativeHypothesis cs + guard $ sided == sidedFromHA ha + +checkTest :: Constraint ComponentSet +checkTest = makeConstraint "check.test" $ \cs -> do + relevance $ guard $ isDerived TestChoice cs + if isInitial TestChoice cs then do + -- initial testType overrules the validTests function + x <- get TestChoice (initials cs) + y <- get TestChoice (derived cs) + guard $ x == y + else do + test <- getTestType TestChoice cs + guard $ test `elem` validTests cs + +checkTestValue, checkTestFormula :: Constraint ComponentSet +checkTestFormula = checkTestValueFor "test-formula" TestFormula +checkTestValue = checkTestValueFor "test-value" TestValue + +checkTestValueFor :: String -> ComponentId -> Constraint ComponentSet +checkTestValueFor s cid = makeConstraint ("check" # s) $ \cs -> do + relevance $ guard $ isDerived cid cs + rel <- getRelation cid cs + let testTypes = inferTestChoices cs + testType <- + case leftHandSide rel of + Var "z" | ZTest `elem` testTypes -> return ZTest + Var "r" | RPearson `elem` testTypes -> return RPearson + Var "F" | Anova `elem` testTypes -> return Anova + Var "chisq" | ChiSquared `elem` testTypes -> return ChiSquared + Var "t" -> + case filter isTTest testTypes of + [] -> fail "test mismatch" + hd:_ -> return hd + _ -> fail "test mismatch" + let tv = rightHandSide rel + sub = getSubstitution cs + tv' <- matchM doubleView (sub |-> tv) + -- special case: when checking test value, first check wether an initial + -- test value is present + case getRhsExpr TestValue (initials cs) of + Just initialTv | cid == TestValue -> do + val <- matchM doubleView initialTv + testDoubles tv' val + _ -> do + target <- case getRhsExpr TestFormula (initials cs) of + Just rhs -> return rhs + _ | testType == ChiSquared -> rightHandSide <$> chiSquaredTestValue cs + _ -> testFormulaFromTest testType (chooseTypeOfTest cs) + target' <- matchM doubleView (sub |-> target) + testDoubles tv' target' + +checkDf :: Constraint ComponentSet +checkDf = makeConstraint "check.df" $ \cs -> do + relevance $ guard $ isDerived Df cs + let cs' = substitute cs + df <- matchM doubleView =<< getExpr Df cs' + case getExpr Df (initials cs) of + Just initialDf -> do + val <- matchM doubleView initialDf + guard $ equalDouble df val + _ -> do + test <- inferTestChoice cs + target <- case test of + ChiSquared -> toExpr <$> chiSquaredDf cs + _ -> degreesOfFreedomFromTest test (chooseTypeOfTest cs) + let sub = getSubstitution cs + target' <- matchM doubleView (sub |-> target) + guard $ equalDouble df target' + +checkDfAnova :: Constraint ComponentSet +checkDfAnova = makeConstraint "check.df-anova" $ \cs -> do + relevance $ guard $ isDerived DfBetween cs || isDerived DfWithin cs + let cs' = substitute cs + (between, within) <- inferDfBetweenWithin cs' + -- check df-between + dfB <- getExpr DfBetween (derived cs) >>= matchM doubleView + unless (equalDouble dfB between) $ fail "df-between" + -- check df-within + dfW <- getExpr DfWithin (derived cs) >>= matchM doubleView + unless (equalDouble dfW within) $ fail "df-within" + +checkGroups :: Constraint ComponentSet +checkGroups = makeConstraint "check.groups" $ \cs -> do + relevance $ guard $ isDerived Groups cs + groups <- getExpr Groups (derived cs) >>= matchM doubleView + -- for now, only support for two groups + guard (groups == 2) + +checkCritical :: Constraint ComponentSet +checkCritical = makeConstraint "check.critical" $ \cs -> do + relevance $ guard $ isDerived Critical cs + unless (derived cs `contains` AlternativeHypothesis) $ + fail "alternative hypothesis missing" + let cs' = substitute cs + critical <- getRelation Critical cs + criticalValue <- matchM doubleView (rightHandSide critical) + let testTypes = inferTestChoices cs + testType <- + case leftHandSide critical of + Var "zcrit" | ZTest `elem` testTypes -> return ZTest + Var "rcrit" | RPearson `elem` testTypes -> return RPearson + Var "Fcrit" | Anova `elem` testTypes -> return Anova + Var "chicrit" | ChiSquared `elem` testTypes -> return ChiSquared + Var "tcrit" -> + case filter isTTest testTypes of + [] -> fail "test mismatch" + hd:_ -> return hd + this -> fail $ "test mismatch " ++ show this + sided <- inferSidedness cs + alpha <- matchM doubleView =<< getExpr SignificanceLevel cs + value' <- + case testType of + Anova -> do + (dfBetween, dfWithin) <- inferDfBetweenWithin cs + computeCriticalF dfBetween dfWithin alpha + ChiSquared -> do + df <- chiSquaredDf cs + computeCriticalChi sided alpha (fromIntegral df) + _ -> do + -- TO DO: use computeCriticalF for Anova + let cs'' = case inferDf cs' of + Just df -> substitute (append Df (CExpr df) cs') + Nothing -> cs' + let df = matchM doubleView =<< getExpr Df cs'' + computeCritical testType sided alpha df + testDoubles value' criticalValue + +checkRejectionCritical :: Constraint ComponentSet +checkRejectionCritical = makeConstraint "check.rejectioncritical" $ \cs -> do + relevance $ guard $ isDerived RejectionCritical cs + unless (derived cs `contains` AlternativeHypothesis) $ + fail "alternative hypothesis missing" + crit <- getRelation RejectionCritical cs + sided <- inferSidedness cs + + let testTypes = inferTestChoices cs + (testType, rel) <- + case withoutAbs (leftHandSide crit) of + Var "z" | ZTest `elem` testTypes -> + return (ZTest, sidedRelation sided (Var "z") (Var "zcrit")) + Var "r" | RPearson `elem` testTypes -> + return (RPearson, sidedRelation sided (Var "r") (Var "rcrit")) + Var "F" | Anova `elem` testTypes -> + return (Anova, sidedRelation sided (Var "F") (Var "Fcrit")) + expr | expr == chiSquared && ChiSquared `elem` testTypes -> + return (ChiSquared, sidedRelation sided chiSquared (Var "chicrit")) + Var "t" -> + case filter isTTest testTypes of + [] -> fail "test mismatch" + hd:_ -> return (hd, sidedRelation sided (Var "t") (Var "tcrit")) + expr -> fail $ "test mismatch" ++ show expr + + unless (withoutAbs (leftHandSide rel) == withoutAbs (leftHandSide crit) + && rightHandSide rel == rightHandSide crit) $ + fail "test mismatch" + unless (leftHandSide rel == leftHandSide crit + && relationType rel == relationType crit) $ + if testType == Anova + then fail "anova mismatch" + else if testType == ChiSquared + then fail "chi-squared mismatch" + else fail "sidedness mismatch" + +withoutAbs :: Expr -> Expr +withoutAbs (Sym s [a]) | s == newSymbol absSymbol = a +withoutAbs expr = expr + +checkPValue :: Constraint ComponentSet +checkPValue = makeConstraint "check.p-value" $ \cs -> do + relevance $ guard $ isDerived PValue cs + unless (derived cs `contains` AlternativeHypothesis) $ + fail "alternative hypothesis missing" + let cs' = substitute cs + t <- matchM doubleView =<< fmap rightHandSide (getRelation TestValue cs') + test <- inferTestChoice cs + pvalue <- matchM doubleView =<< getExpr PValue cs + when (pvalue < 0 || pvalue > 1) $ + fail "value not a probability" + sided <- inferSidedness cs + let cs'' = case inferDf cs' of + Just df -> substitute (append Df (CExpr df) cs') + Nothing -> cs' + let df = matchM doubleView =<< getExpr Df cs'' + pvalue' <- computePValue test sided t df + -- compare with 3 decimals + testDoublesWith 0.00055 pvalue' pvalue $ + if derived cs `contains` TestValue + then empty + else fail "TestValue missing" + +checkConclusionCritical :: Constraint ComponentSet +checkConclusionCritical = makeConstraint "check.conclusion-critical" $ \cs -> do + relevance $ guard $ isDerived ConclusionCritical cs + rej <- inferRejectionCritical cs + -- to do: rejection critical is added to component set only to get the substituted relation + let cs' = substitute (append RejectionCritical (CRelation rej) cs) + rejection <- getRelation RejectionCritical cs' + lhs <- matchM doubleView $ leftHandSide rejection + rhs <- matchM doubleView $ rightHandSide rejection + let result = eval (relationType rejection) lhs rhs + conclusion <- getConclusion ConclusionCritical cs + guard $ result == conclusion + +checkConclusionHypotheses :: Constraint ComponentSet +checkConclusionHypotheses = makeConstraint "check.conclusion-hypotheses" $ \cs -> do + relevance $ guard $ isDerived ConclusionHypotheses cs + concl <- getRejectionHypotheses ConclusionHypotheses cs + b1 <- case fmap relationType (inferConclusionPValue cs) of + Just LessThanOrEqualTo + | concl == RejectH0 -> return True + Just GreaterThan + | concl == DontRejectH0 -> return True + Just _ -> fail "conclusion mismatch pvalue" + Nothing -> return False + b2 <- case inferConclusionCritical cs of + Just False + | concl == DontRejectH0 -> return True + | concl `elem` [AcceptH0, RejectH1] -> fail "convention DontRejectH0" + Just True + | concl == RejectH0 -> return True + | concl `elem` [AcceptH1, DontRejectH1] -> fail "convention RejectH0" + Just _ -> fail "conclusion mismatch critical" + Nothing -> return False + guard (b1 || b2) + +checkStandardError :: Constraint ComponentSet +checkStandardError = makeConstraint "check.standard-error" $ \cs -> do + relevance $ guard $ isDerived StandardError cs + rel <- getRelation StandardError cs + n <- matchM doubleView =<< getExpr SampleSize cs + case (leftHandSide rel, getRhsExpr PopulationSdev cs, getRhsExpr SampleSdev cs) of + (Var "sigmaM", Just expr, Nothing) -> do + psdev <- matchM doubleView expr + se <- matchM doubleView (rightHandSide rel) + testDoubles se (psdev/sqrt n) + (Var "SEM", Nothing, Just expr) -> do + sdev <- matchM doubleView expr + se <- matchM doubleView (rightHandSide rel) + testDoubles se (sdev/sqrt n) + _ -> fail "standard error mismatch" + +checkObservedTotals :: Constraint ComponentSet +checkObservedTotals = makeConstraint "check.observed-totals" $ \cs -> do + relevance $ guard $ + all (derived cs `contains`) [ObservedRowTotals, ObservedColumnTotals, ObservedTotal] + rowTotals <- getExpr ObservedRowTotals cs >>= fromExpr + columnTotals <- getExpr ObservedColumnTotals cs >>= fromExpr + total <- getExpr ObservedTotal cs >>= fromExpr + table <- getTable ObservedFrequencies cs + let (expRowTotals, expColumnTotals, expTotal) = computeTotals table + unless (rowTotals == expRowTotals) $ fail "rows" + unless (columnTotals == expColumnTotals) $ fail "columns" + unless (total == expTotal) $ fail "total" + +checkExpectedFrequencies :: Constraint ComponentSet +checkExpectedFrequencies = makeConstraint "check.expected-frequencies" $ \cs -> do + relevance $ guard $ derived cs `contains` ExpectedFrequencies + expectedFrequencies <- getExpr ExpectedFrequencies cs >>= matchDoubleTable + observed <- getTable ObservedFrequencies cs + let totals = computeTotals observed + computedFrequencies = computeExpectedFrequencies totals + guard (equalTableDouble expectedFrequencies computedFrequencies) + +matchDoubleTable :: MonadPlus m => Expr -> m [[Double]] +matchDoubleTable = fromExpr >=> mapM (mapM (matchM doubleView)) + +equalTableDouble :: [[Double]] -> [[Double]] -> Bool +equalTableDouble = equalListBy (equalListBy equalDouble) + +equalListBy :: (a -> a -> Bool) -> [a] -> [a] -> Bool +equalListBy eq xs ys = length xs == length ys && and (zipWith eq xs ys) + +---------------------------------------------------------- + +testDoubles :: Double -> Double -> Result () +testDoubles value target = testDoublesWith defaultDelta value target empty + +testDoublesWith :: Double -> Double -> Double -> Result () -> Result () +testDoublesWith delta value target resultNotEqual = + unless (equalDoubleWith delta value target) $ + if equalDoubleWith (delta*10) value target + then fail "almost equal" + else resultNotEqual + +equalDouble :: Double -> Double -> Bool +equalDouble = equalDoubleWith defaultDelta + +equalDoubleWith :: Double -> Double -> Double -> Bool +equalDoubleWith delta x y = abs (x - y) < delta + +defaultDelta :: Double +defaultDelta = 0.0055 + +-- afronden op 3 cijfers achter de komma +-- rounded :: Double -> Double +-- rounded x = fromInteger (round (x * 1000)) / 1000.0
+ src/Domain/Hypothesis/Examples.hs view
@@ -0,0 +1,450 @@+-----------------------------------------------------------------------------+-- Copyright 2020, Ideas project team. This file is distributed under the+-- terms of the Apache License 2.0. For more information, see the files+-- "LICENSE.txt" and "NOTICE.txt", which are included in the distribution.+-----------------------------------------------------------------------------+module Domain.Hypothesis.Examples+ {- ( hypothesisExamples+ , ex_dwo2+ , ex_dwo3+ , ex_simple+ ) -} where++import Domain.Math.Data.Relation+import Domain.Math.Expr+import Domain.Statistics.ComponentSet+import Ideas.Common.Exercise+import Prelude hiding ((^))++hypothesisExamples :: Examples ComponentSet+hypothesisExamples = examplesWithDifficulty $+ [ (Easy, ex_simple)+ , (Easy, ex_omi1)+ , (Easy, ex_omi2)+ , (Easy, ex_omi3)+ , (Medium, ex_econ)+ , (Medium, ex_dwo1)+ , (Medium, ex_dwo2)+ , (Medium, ex_dwo3)+ , (Medium, ex_dwo4)+ , (Medium, ex_dwo5)+ , (Medium, ex_dwo6)+ , (Medium, ex_dwo7)+ , (Medium, sietske_1)+ , (Medium, sietske_2)+ , (Medium, sietske_3)+ , (Difficult, ex_bio')+ ] +++ [ (Medium, cs) | cs <- [mayPilotA, mayPilotB, mayPilotC, mayPilotD] ]++cExprDouble :: Double -> Component+cExprDouble = CExpr . fromDouble++--------------------------------------------------------------------------------+-- * Pilot exercises (May 2017) * --+--------------------------------------------------------------------------------++{- Montarello and Martins (2005) discovered that children from the 7th group+ were able to solve more difficult maths exercises if very simple exercises+ were added to the exam. To research this phenomenon a researcher picks a+ standardised exam of which the scores have a gaussian distribution with++ mu = 50,+ sigma = 8,++ and adds very easy exercises to it. They give the exam to 40 children. The+ mean outcome of the test is M = 52. Is this result enough to conclude that+ adding the easy exercises improves the childrens prestation? Assume a+ significance level of alpha = 0.05 -}++-- Note: Claim assumed that `improvement' means mu_after - mu_before > 0+mayPilotA :: ComponentSet+mayPilotA = initialSet+ [ -- ( Claim , CRelation $ Var "mu" .>. fromDouble 50 )+ ( AlternativeHypothesis, CRelation $ Var "mu" .>. fromDouble 50)+ , ( SampleMean , cExprDouble 52 )+ , ( PopulationMean , cExprDouble 50 ) -- I think that different means will+ , ( PopulationSdev , cExprDouble 8 ) -- cause all sorts of trouble with+ , ( SampleSize , CExpr 40 ) -- the strategy.+ , ( SignificanceLevel , cExprDouble 0.05 ) ]++{- It has been showed by Evans, Pelham, Smith et al. (2001) that the Ritalin+ medicin improves the attention span of children with ADHD and that their+ school results improve.++ A researcher draws a sample of N = 16 children that have been diagnosed with+ ADHD to demonstrate that the medicine truly works. Their attention span is+ measured before and after the medicine has been used. The results show an+ average improvement of the attention span of M = 3.13 minutes with a standard+ deviation of s = 5.39. Are the results of this experiment enough to conclude+ that the the medicine improves the attention span of children with ADHD?+ Assume a significance level of alpha = 0.01. -}++-- Note: mu_d means mu_after - mu_before+mayPilotB :: ComponentSet+mayPilotB = initialSet+ [ -- ( Claim , CRelation $ Var "mu" .>. fromDouble 0)+ ( AlternativeHypothesis, CRelation $ Var "mu" .>. fromDouble 0)+ , ( SampleSize , CExpr 16 )+ , ( SampleMean , cExprDouble 3.13 )+ , ( SampleSdev , cExprDouble 5.39 )+ , ( SignificanceLevel , cExprDouble 0.01 ) + , ( TestChoice , CChoice (TestType TTestPaired))+ ]++{- A researcher compares two treatments during a design with independent groups.+ Group A has 18 participants and group V has 16 participants. The means for the+ groups are found to be M = 12.14 for group A and M = 16.17 for group B. The+ pooled variance s^2 is approximated to be 7.74. The researcher wants to know+ if the treatments will leaad to a different result. Assume a significance+ level of alpha = 0.05. -}++mayPilotC :: ComponentSet+mayPilotC = initialSet+ [ -- ( Claim , CRelation $ Var "mu1" ./=. Var "mu2")+ ( AlternativeHypothesis, CRelation $ Var "mu1" ./=. Var "mu2")+ , ( One SampleSize , CExpr 18 )+ , ( Two SampleSize , CExpr 16 )+ , ( One SampleMean , cExprDouble 12.14 )+ , ( Two SampleMean , cExprDouble 16.70 )+ , ( PooledVariance , cExprDouble 7.74 )+ , ( SignificanceLevel , cExprDouble 0.05 ) + ]++{- The mean annual income for a randomised sample of 51 househoulds turns out to+ be 23 613 euros. The standard deviation of the sample is 658 euros. You may+ assume that annual income is normally distributed. Research the claim that+ mean annual income in the population is lower than 24000. Assume a+ significance level of alpha = 0.10. -}++mayPilotD :: ComponentSet+mayPilotD = initialSet+ [ -- ( Claim , CRelation $ Var "mu" .<. fromDouble 24000)+ ( AlternativeHypothesis, CRelation $ Var "mu" .<. fromDouble 24000)+ , ( SampleSize , CExpr 51 )+ , ( SampleMean , cExprDouble 23613 )+ , ( PopulationMean , cExprDouble 24000 )+ , ( SampleSdev , cExprDouble 658 )+ , ( SignificanceLevel , cExprDouble 0.10 ) + ]++-- Example from an exam+ex_econ :: ComponentSet+ex_econ = initialSet+ [ (SampleSize, CExpr 61)+ , (SampleMean, CExpr $ fromDouble 21023.0)+ , (SampleSdev, CExpr $ fromDouble 634.0)+ , (AlternativeHypothesis, CRelation $ Var "mu" .<=. 22000)+ ]++-- Example for an hypothesis test+ex_simple :: ComponentSet+ex_simple = initialSet [+ (AlternativeHypothesis, CRelation $ Var "mu" .>=. fromDouble 98.6),+ (SampleSize, CExpr 101),+ (SampleMean, CExpr $ fromDouble 98.9),+ (SampleSdev, CExpr $ fromDouble 0.6)]++-- Examples from the OMI course+ex_omi1 :: ComponentSet+ex_omi1 = initialSet [+ (AlternativeHypothesis, CRelation $ Var "mu" ./=. fromDouble 20.0),+ (SampleSize, CExpr 10),+ (SampleMean, CExpr $ fromDouble 22.3),+ (SampleSdev, CExpr $ fromDouble 3.65)+ ]++ex_omi2 :: ComponentSet+ex_omi2 = initialSet [+ (AlternativeHypothesis, CRelation $ Var "mu" .>=. fromDouble 20.0),+ (SampleSize, CExpr 10),+ (SampleMean, CExpr $ fromDouble 22.3),+ (SampleSdev, CExpr $ fromDouble 3.65)+ ]++ex_omi3 :: ComponentSet+ex_omi3 = initialSet [+ (AlternativeHypothesis, CRelation $ Var "mu1" .==. Var "mu2"),+ (One SampleSize, CExpr 46),+ (Two SampleSize, CExpr 56),+ (One SampleMean, CExpr $ fromDouble 3.8698),+ (Two SampleMean, CExpr $ fromDouble 4.5819),+ (One SampleSdev, CExpr $ fromDouble 1.6714),+ (Two SampleSdev, CExpr $ fromDouble 1.5216),+ (PooledVariance, CExpr $ fromDouble 2.531)+ ]++-- Example from an exam simplified+ex_bio' :: ComponentSet+ex_bio' = initialSet [+ (AlternativeHypothesis, CRelation $ Var "mu1" .==. Var "mu2"),+ (One SampleSize, CExpr 9),+ (Two SampleSize, CExpr 9),+ (One SampleMean, CExpr $ fromDouble 19.189),+ (Two SampleMean, CExpr $ fromDouble 28.067),+ (One SampleSdev, CExpr $ fromDouble 6.36),+ (Two SampleSdev, CExpr $ fromDouble 6.12),+ (PooledVariance, CExpr $ fromDouble 38.843)]++-- Examples from the DWO+ex_dwo1 :: ComponentSet+ex_dwo1 = ex_econ++ex_dwo2 :: ComponentSet+ex_dwo2 = initialSet [+ (SampleSize, CExpr 68),+ (SampleMean, CExpr $ fromDouble 23035.0),+ (SampleSdev, CExpr $ fromDouble 658.0),+ (AlternativeHypothesis, CRelation $ Var "mu" .>=. 24000),+ (SignificanceLevel, CExpr $ fromDouble 0.10)]++ex_dwo3 :: ComponentSet+ex_dwo3 = initialSet [+ (SampleSize, CExpr 61),+ (SampleMean, CExpr $ fromDouble 21023.0),+ (SampleSdev, CExpr $ fromDouble 634.0),+ (AlternativeHypothesis, CRelation $ Var "p" ./=. fromDouble 0.25),+ (Proportion, CExpr $ fromDouble 0.21)+-- fix me (TestStatistic "t", CExpr $ (Var "p" - Var "p0") / sqrt ((Var "p0") * (1.0 - Var "p0") / (Var "n")))+ ]++ex_dwo4 :: ComponentSet+ex_dwo4 = initialSet [+ (SampleSize, CExpr 72),+ (SampleMean, CExpr $ fromDouble 24061.0),+ (SampleSdev, CExpr $ fromDouble 663.0),+ (AlternativeHypothesis, CRelation $ Var "p" ./=. fromDouble 0.3),+ (Proportion, CExpr $ fromDouble 0.22),+ --(TestStatistic "t", CExpr $ (Var "p" - Var "p0") / sqrt ((Var "p0") * (1.0 - Var "p0") / (Var "n"))),+ (SignificanceLevel, CExpr $ fromDouble 0.01)]++ex_dwo5 :: ComponentSet+ex_dwo5 = initialSet [+ (SampleSize, CExpr 100),+ (SampleMean, CExpr $ fromDouble 24061.0),+ (SampleSdev, CExpr $ fromDouble 663.0),+ (AlternativeHypothesis, CRelation $ Var "p1" .>=. Var "p2"),+ (One Proportion, CExpr $ fromDouble 0.3),+ (Two Proportion, CExpr $ fromDouble 0.1),+ (Other "p0", CExpr $ (Var "p1" + Var "p2") / 2),+ (Other "d0", CExpr $ fromDouble 0.0)+-- fix me (TestStatistic "t", CExpr $ (Var "p1" - Var "p2" - Var "d0") / (sqrt $ (Var "p0") * (1.0 - Var "p0") / ((Var "n") / 2)))+ ]++ex_dwo6 :: ComponentSet+ex_dwo6 = initialSet [+ (SampleSize, CExpr 100),+ (SampleMean, CExpr $ fromDouble 24061.0),+ (SampleSdev, CExpr $ fromDouble 663.0),+ (AlternativeHypothesis, CRelation $ Var "p1" ./=. Var "p2"),+ (One Proportion, CExpr $ fromDouble 0.3),+ (Two Proportion, CExpr $ fromDouble 0.1),+ (Other "p0", CExpr $ (Var "p1" + Var "p2") / 2),+ (Other "d0", CExpr $ fromDouble 0.0)+-- fix me (TestStatistic "t", CExpr $ (Var "p1" - Var "p2" - Var "d0") / (sqrt $ (Var "p0") * (1.0 - Var "p0") / ((Var "n") / 2)))+ ]++ex_dwo7 :: ComponentSet+ex_dwo7 = initialSet [+ (SampleSize, CExpr 490),+ (SampleMean, CExpr $ fromDouble 83.0),+ (Correlation, CExpr $ fromDouble (-0.557)),+ (Df, CExpr $ Var "n" - 2),+ (AlternativeHypothesis, CRelation $ Var "r" .>=. fromDouble 0.0)+ ]+++{- "Average income" example+ -+ - Steekproefgrootte n = 61+ - Steekproefgemiddelde m = 21023+ - Steekproefstandaardafwijking s = 634+ - “Je mag aannemen dat inkomen normaal is verdeeld”+ - Significantieniveau α = 0,05.+ - Claim: Het populatiegemiddelde μ is kleiner dan 22000.+ -}++sietske_1 :: ComponentSet+sietske_1 = initialSet [+ (SampleSize, CExpr 61),+ (SampleMean, CExpr 21023),+ (SampleSdev, CExpr 634),+ (SignificanceLevel, CExpr $ fromDouble 0.05),+ (AlternativeHypothesis, CRelation $ Var "mu" .<=. 22000)+ ]++{- "Car factory" example+ -+ - Steekproefgrootte: n_1 = 100, n_2 = 100+ - Steekproefproportie 1: p_A = 0.3+ - Steekproefproportie 2: p_B = 0.1+ - Significantieniveau: alpha = 0.05.+ - p_0 = (p_A + p_B) / 2+ - Verschil tussen populatieproporties onder nulhypothese D_0 = 0+ - Toetsingsgrootheid+ - Claim: Populatieproportie P_A is groter dan populatieproportie P_B+ -}+sietske_2 :: ComponentSet+sietske_2 = initialSet+ [ (SampleSize , CExpr 100)+ , (One Proportion , CExpr $ fromDouble 0.3)+ , (Two Proportion , CExpr $ fromDouble 0.1)+ , (SignificanceLevel, CExpr $ fromDouble 0.05)+ , (Other "p0" , CExpr $ (Var "p1" + Var "p2") / 2)+ , (Other "d0" , CExpr $ fromDouble 0.0)+-- fix me+-- , (TestStatistic "t",+-- CExpr $ (Var "p1" - Var "p2" - Var "d0") / (sqrt $ (Var "p0") *+-- (1.0 - Var "p0") / ((Var "n") / 2)))+ , (AlternativeHypothesis , CRelation $ Var "p1" .>=. Var "p2")+ ]++{- "Shop owner" example+ -+ - Steekproefgrootte n_x=20 en n_y=20+ - Steekproefgemiddelde 1: = 85+ - Steekproefgemiddelde 2: = 63+ - Steekproefstandaardafwijking 1: s_x = 11+ - Steekproefstandaardafwijking 2: s_y = 11+ - Aanname dat normaal verdeeld+ - Verschil tussen populatiegrootheden onder nulhypothese D_0=0+ - Toetsingsgrootheid+ - Significantieniveau α = 0,05.+ - Aantal vrijheidsgraden DF=n_x+n_y-2+ - Claim: Populatiegemiddelde mu_A is ongelijk aan populatiegemiddelde mu_B.+ -}+sietske_3 :: ComponentSet+sietske_3 = initialSet [+ (SampleSize, CExpr 100),+ (One SampleMean, CExpr 85),+ (Two SampleMean, CExpr 63),+ (One SampleSdev, CExpr 11),+ (Two SampleSdev, CExpr 11),+-- fix me(TestStatistic "t", CExpr $ (mean1 - mean2 - (Var "d0")) / (sqrt((sd1 ^ 2 / (Var "n")) + (sd2 ^ 2 / (Var "n")))) ),+ (SignificanceLevel, CExpr $ fromDouble 0.05),+ -- (Df, CExpr $ (Var "n1) + (Var "n2") - 2),+ (AlternativeHypothesis, CRelation $ Var "muA" ./=. Var "muB")+ ]+ {-+ where+ mean1 = toExpr (One SampleMean)+ mean2 = toExpr (Two SampleMean)+ sd1 = toExpr (One SampleSdev)+ sd2 = toExpr (Two SampleSdev) -}++--------------------------------------------------------------------------------+-- Opgaven Sociale Wetenschappen++opgavenSW :: [ComponentSet]+opgavenSW = [ -- december 2017+ opgave_3_4, opgave_3_6, opgave_4_10, opgave_4_11, opgave_5_3, opgave_5_6+ -- februari 2018+ , opgave_1_5, opgave_3_5, opgave_5_4+ ]++--------------------------------------------------------------------------------+-- Opgaven Sociale Wetenschappen (december 2017)++opgave_3_4 :: ComponentSet+opgave_3_4 = initialSet + [ (AlternativeHypothesis, CRelation $ Var "mu" .>. 100)+ , (PopulationSdev, CRelation $ Var "sigma" .==. 18)+ , (SampleMean, CExpr 104)+ , (SampleSize, CExpr 36)+ , (SignificanceLevel, CExpr 0.01)+ ]++opgave_3_6 :: ComponentSet+opgave_3_6 = initialSet + [ (AlternativeHypothesis, CRelation $ Var "mu" .<. 4.9)+ , (PopulationSdev, CRelation $ Var "sigma" .==. 0.84)+ , (SampleMean, CExpr 4.4)+ , (SampleSize, CExpr 16)+ , (SignificanceLevel, CExpr 0.05)+ ]++opgave_4_10 :: ComponentSet+opgave_4_10 = initialSet + [ (AlternativeHypothesis, CRelation $ Var "mu" ./=. 0)+ , (SampleSdev, CRelation $ Var "s" .==. 1.50)+ , (SampleMean, CExpr 1.28)+ , (SampleSize, CExpr 25)+ , (TestChoice, CChoice $ TestType TTestPaired)+ , (SignificanceLevel, CExpr 0.05)+ ]++opgave_4_11 :: ComponentSet+opgave_4_11 = initialSet + [ (AlternativeHypothesis, CRelation $ Var "mu" ./=. 0)+ , (SampleSdev, CRelation $ Var "s" .==. 2.45)+ , (SampleMean, CExpr 2)+ , (SampleSize, CExpr 6)+ , (TestChoice, CChoice $ TestType TTestPaired)+ , (SignificanceLevel, CExpr 0.05)+ ]++opgave_5_3 :: ComponentSet+opgave_5_3 = initialSet+ [ (AlternativeHypothesis, CRelation $ Var "mu1" ./=. Var "mu2")+ , (Df, CExpr 30)+ , (TestValue, CRelation $ Var "t" .==. 2.085)+ , (TestChoice, CChoice $ TestType TTestTwo)+ , (SignificanceLevel, CExpr 0.05)+ ]++opgave_5_6 :: ComponentSet+opgave_5_6 = initialSet+ [ (AlternativeHypothesis, CRelation $ Var "mu1" ./=. Var "mu2")+ , (Df, CExpr 119.50)+ , (TestValue, CRelation $ Var "t" .==. 3.379)+ , (TestChoice, CChoice $ TestType TTestTwo)+ , (SignificanceLevel, CExpr 0.05)+ ]+ +--------------------------------------------------------------------------------+-- Opgaven Sociale Wetenschappen (februari 2018)++opgave_1_5 :: ComponentSet+opgave_1_5 = initialSet+ [ (AlternativeHypothesis, CRelation $ Var "rho" ./=. 0)+ , (SignificanceLevel, CExpr 0.01)+ , (TestValue, CRelation $ Var "r" .==. 0.835)+ , (SampleSize, CExpr 9) + , (TestChoice, CChoice $ TestType RPearson)+ ]++opgave_3_5 :: ComponentSet+opgave_3_5 = initialSet+ [ (AlternativeHypothesis, CRelation $ Var "mu1" ./=. Var "mu2")+ , (TestValue, CRelation $ Var "F" .==. 4.00)+ , (SignificanceLevel, CExpr 0.05)+ , (SampleSize, CExpr 40)+ , (TestChoice, CChoice $ TestType Anova)+ ]+ +opgave_5_4 :: ComponentSet+opgave_5_4 = initialSet+ [ (SignificanceLevel, CExpr 0.01)+ , (TestChoice, CChoice $ TestType ChiSquared)+ , (ObservedFrequencies, CExpr $ toExpr [[18 :: Int, 4, 2], [4, 17, 15]])+ ]++--------------------------------------------------------------------------------+-- Opgaven Economie (februari 2018)++opgavenEconomie :: [ComponentSet]+opgavenEconomie = [opgave_4_18]+ +opgave_4_18 :: ComponentSet+opgave_4_18 = initialSet+ [ (AlternativeHypothesis, CRelation $ Var "mu1" ./=. Var "mu2")+ , (SignificanceLevel, CExpr 0.05)+ , (One SampleSize, CExpr 20) + , (Two SampleSize, CExpr 20)+ , (One SampleMean, CExpr 85)+ , (Two SampleMean, CExpr 63)+ , (One SampleSdev, CExpr 11)+ , (Two SampleSdev, CExpr 11)+ , (TestFormula, CRelation $ Var "t" .==. + (Var "samplemean1" - Var "samplemean2") / sqrt ((Var "samplesd1" ^ 2) / Var "n1" + (Var "samplesd2" ^ 2) / Var "n2"))+ ]
+ src/Domain/Hypothesis/Exercises.hs view
@@ -0,0 +1,178 @@+-----------------------------------------------------------------------------+-- Copyright 2020, Ideas project team. This file is distributed under the+-- terms of the Apache License 2.0. For more information, see the files+-- "LICENSE.txt" and "NOTICE.txt", which are included in the distribution.+-----------------------------------------------------------------------------+module Domain.Hypothesis.Exercises (hypothesisExercise) where++import Control.Monad+import qualified Data.Map as M+import Data.List+import Data.Maybe+import Domain.Math.Data.Relation+import Domain.Hypothesis.Examples+import Domain.Hypothesis.Constraints+import Domain.Hypothesis.BuggyRules+import Domain.Hypothesis.Rules+import Domain.Hypothesis.Strategies+import Domain.Statistics.ComponentSet+import Domain.Statistics.Parser+import Ideas.Common.Library hiding (Predicate)+import Domain.Math.Expr.Data+import Ideas.Encoding.OpenMathSupport+import Ideas.Utils.Uniplate+import Ideas.Text.XML hiding (children)+import Prelude hiding (until)++-- Debug functions+css :: [ComponentSet]+css = mapMaybe fromContext $ concatMap terms $+ mapMaybe (defaultDerivation hypothesisExercise) (examplesAsList hypothesisExercise)++_testPP :: IO ()+_testPP = forM_ css $ \x ->+ let s = prettyPrinter hypothesisExercise x in+ case parser hypothesisExercise s of+ Right y | x == y -> putChar '.'+ | otherwise -> error $ show (x, y)+ Left msg -> error $ s ++ "\n" ++ msg++_testOM :: IO ()+_testOM = forM_ css $ \x ->+ case toOpenMath hypothesisExercise x of+ Right omobj ->+ case fromOpenMath hypothesisExercise omobj of+ Just y | x == y -> putChar '.'+ | otherwise -> error $ "Not the same:\n" ++ show x ++ "\n" ++ show y+ Nothing -> error $ show omobj+ Left msg -> error $ show x ++ "\n" ++ msg++_see :: Int -> IO ()+_see n = printDerivation hypothesisExercise (examplesAsList hypothesisExercise !! n)++_testje :: IO ()+_testje = printDerivation hypothesisExercise opgave_5_4++_save :: IO ()+_save = do+ let pilot = [("A", mayPilotA), ("B", mayPilotB), ("C", mayPilotC), ("D", mayPilotD)]+ ex = hypothesisExercise+ forM_ pilot $ \(n, cs) -> do+ writeFile ("Task" ++ n ++ ".txt") $ showDerivation ex cs+ writeFile ("Task" ++ n ++ "-OM.txt") $ showDerivationOM ex cs++showDerivationOM :: Exercise a -> a -> String+showDerivationOM ex a =+ case defaultDerivation ex a of+ Just d -> show $ biMap fst f d+ Nothing -> "no derivation"+ where+ errorXML = makeXML "error" mempty+ f ctx = fromMaybe errorXML $ do+ x <- fromContext ctx+ omobj <- toOpenMath ex x+ return (toXML omobj)++------------------------------++hypothesisExercise :: Exercise ComponentSet+hypothesisExercise = emptyExercise+ { exerciseId = describe "Hypothesis testing" $ newId "hypothesis"+ , prettyPrinter = show+ , status = Experimental+ , parser = parseComponentSet+ , strategy = liftToContext hypothesisStrategy+ , extraRules = fmap liftToContext buggyRules+ , ruleOrdering = ruleOrderingWith highPriorityRules+ , examples = examplesWithDifficulty [ (Medium, cs) | cs <- opgavenSW ++ opgavenEconomie ] -- hypothesisExamples+ , equivalence = withoutContext eqComponentSet+ , similarity = withoutContext similarComponentSet+ , constraints = map liftToContext hypothesisConstraints+ , hasTypeable = useTypeable+ , ready = predicate (hasConclusionHypotheses <&&> checkConstraints)+ , hasTermView = Just termView+ }+ +highPriorityRules :: [Rule ComponentSet]+highPriorityRules = + [ addHypothesesRule+ , lookupTValueRule, lookupZValueRule+ , computePValueTTest, computePValueZTest+ , addTestValueRule+ , hypothesesConclusionCriticalRule, hypothesesConclusionPValueRule+ -- buggy rules+ , buggyTValueSided, buggyTValueTestValue -- Sietske: higher priority than corresponding buggy alpha rule+ , buggyRValueSided, buggyRValueTestValue -- idem+ , buggyZValueSided, buggyZValueTestValue+ , buggyChiValueSided, buggyChiValueTestValue+ ]++----------------------------------------------------------+-- Equivalence+eqComponentSet :: ComponentSet -> ComponentSet -> Bool+eqComponentSet x y = compareIntials+ && checkConstraints x+ && checkConstraints y+ where+ compareIntials =+ M.fromList (toList (initials x)) == M.fromList (toList (initials y))++checkConstraints :: ComponentSet -> Bool+checkConstraints cs =+ all (`checkConstraint` cs) hypothesisConstraints+ +checkConstraint :: Constraint ComponentSet -> ComponentSet -> Bool+checkConstraint p cs + | not (isRelevant p cs) = True+ | otherwise = isSatisfied p cs++----------------------------------------------------------+-- Similarity++-- To do: not all component ids compare doubles with the same precision+-- (e.g. pvalue uses 3 decimals, default is 2)+similarComponentSet :: ComponentSet -> ComponentSet -> Bool+similarComponentSet cs1 cs2 = sorted cs1 `eqList` sorted cs2+ where+ sorted = map f . sortOn fst . toList+ where+ f (n, CExpr e) = (n, CExpr (normalizeExpr e))+ f (n, CRelation r) = (n, CRelation (fmap normalizeExpr r))+ f (n, c) = (n, c)+ + eqList :: [(ComponentId, Component)] -> [(ComponentId, Component)] -> Bool+ eqList xs ys = length xs == length ys && all (\((a, x), (b, y)) -> a == b && x `eqComponent` y) (zip xs ys)++ eqComponent :: Component -> Component -> Bool+ eqComponent (CExpr x) (CExpr y) = eqExpr x y+ eqComponent (CRelation x) (CRelation y) = + eqExpr (leftHandSide x) (leftHandSide y) && eqExpr (rightHandSide x) (rightHandSide y)+ && relationType x == relationType y+ eqComponent x y = x == y+ + eqExprs :: [Expr] -> [Expr] -> Bool+ eqExprs [] [] = True+ eqExprs (x:xs) (y:ys) = eqExpr x y && eqExprs xs ys+ eqExprs _ _ = False+ + eqExpr :: Expr -> Expr -> Bool+ eqExpr (Number x) (Number y) = equalDouble x y+ eqExpr x y = + case (getFunction x, getFunction y) of + (Just (f, xs), Just (g, ys)) -> f == g && xs `eqExprs` ys+ _ -> x == y+++-- we vouwen hier de definities uit zodat we formules kunnen herkennen+-- afronden op 2 decimalen: beter is om expressies component-gewijs te +-- vergelijken en daarbij de hulpfunctie 'equal' te gebruiken+normalizeExpr :: Expr -> Expr+normalizeExpr (Var "sigmaM") = Var "sigma" / sqrt (Var "n")+normalizeExpr (Var "SEM") = Var "s" / sqrt (Var "n")+normalizeExpr (Nat n) = Number (fromInteger n)+normalizeExpr expr = descend normalizeExpr expr++-- | Test whether the conclusionhypothesis is present in the componentset. If+-- the predicate holds then he exercise is considered /done/.+hasConclusionHypotheses :: ComponentSet -> Bool+hasConclusionHypotheses cs = contains cs ConclusionHypotheses
+ src/Domain/Hypothesis/Rules.hs view
@@ -0,0 +1,666 @@+-----------------------------------------------------------------------------+-- Copyright 2020, Ideas project team. This file is distributed under the+-- terms of the Apache License 2.0. For more information, see the files+-- "LICENSE.txt" and "NOTICE.txt", which are included in the distribution.+-----------------------------------------------------------------------------+module Domain.Hypothesis.Rules+ ( addAlphaRule, addDfRule, addDfBetweenWithin+ , addH0FromHARule, addH0FromHAEqualSignRule, addHARule+ , addHypothesesRule, addHypothesesChiSquaredRule+ , addObservedTotals, addExpectedFrequencies+ , addConclusionPValueRule, addRejectionRule+ , addTestFormulaRule, addTestValueRule+ , chooseTTestPairedRule, chooseTTestRule, chooseTTestTwoRule+ , chooseZTestRule, chooseRPearsonRule, chooseAnovaRule, chooseChiSquaredRule+ , computePValueTTest, computePValueZTest+ , criticalConclusionRule+ , hypothesesConclusionCriticalRule, hypothesesConclusionPValueRule+ , addStandardErrorSigma, addStandardErrorSD+ , determineSided+ , lookupTValueRule, lookupZValueRule, lookupRValueRule, lookupFValueRule+ , lookupChiValueRule+ -------------+ , inferSidedness, inferTestChoice, inferTestChoices, inferRejectionCritical+ , inferConclusionCritical, inferConclusionPValue, inferDf, inferVar+ , inferTestFormula, inferCriticalZWith, inferCriticalTWith+ , inferCriticalRWith, inferCriticalFWith, inferCriticalChiWith+ , inferTestValue, inferDfBetweenWithin, chiSquaredDf, chiSquaredTestValue+ , computeTotals, computeExpectedFrequencies, getTable+ ) where++import Control.Monad+import Data.List+import Data.Maybe+import Domain.Hypothesis.Common+import Domain.Math.Data.Relation+import Domain.Math.Expr hiding ((.*.), (./.), (^))+import Domain.Math.Numeric.Views+import Domain.Statistics.ComponentSet+import Domain.Statistics.Views+import Ideas.Common.Library++----------------------------------------------------------+-- Rules for determining the confidence level++addAlphaRule :: Rule ComponentSet+addAlphaRule =+ describe "Rule for adding the alpha component" .+ makeRule "component.alpha" $ f+ where+ f :: ComponentSet -> Maybe ComponentSet+ f cs = do+ guard (derived cs `doesNotContain` SignificanceLevel)+ let alpha = case get SignificanceLevel (initials cs) of+ Just (CExpr a) -> a+ _ -> toExpr pickAlpha+ return $ append SignificanceLevel (CExpr alpha) cs++----------------------------------------------------------+-- Rules for constructing the hypotheses++determineSided :: Rule ComponentSet+determineSided =+ describe "Rule for determine one-/two-sided testing" .+ makeRule "component.sided" $ f+ where+ f :: ComponentSet -> Maybe ComponentSet+ f cs = do+ guard (cs `doesNotContain` Sidedness)+ sided <- + case getTestType TestChoice (initials cs) of+ Just Anova -> return RightSided+ _ -> do+ ha <- getRelation AlternativeHypothesis cs+ return (sidedFromHA ha)+ return $ append Sidedness (CChoice $ Sided sided) cs++addHypothesesRule :: Rule ComponentSet+addHypothesesRule = + describe "Add null hypothesis and alternative hypothesis, in one step" $+ makeRule "component.hypotheses" f+ where+ f :: ComponentSet -> [ComponentSet]+ f = applyAll $ (addH0FromHARule ./. addH0FromHAEqualSignRule) .*. addHARule++addH0FromHARule :: Rule ComponentSet+addH0FromHARule =+ describe "Rule for adding the H0 component based on HA" .+ makeRule "component.h0-from-ha" $ f+ where+ f :: ComponentSet -> Maybe ComponentSet+ f cs = do+ guard (cs `doesNotContain` NullHypothesis)+ ha <- getRelation AlternativeHypothesis cs+ return . append NullHypothesis+ (CRelation $ h0FromHA ha) $ cs++addH0FromHAEqualSignRule :: Rule ComponentSet+addH0FromHAEqualSignRule =+ describe "Rule for adding the H0 component based on HA; use equal sign (by convention)" .+ makeRule "component.h0-from-ha-eq" $ f+ where+ f :: ComponentSet -> Maybe ComponentSet+ f cs = do+ guard (cs `doesNotContain` NullHypothesis)+ ha <- getRelation AlternativeHypothesis cs + return . append NullHypothesis+ (CRelation $ h0FromHAEqualSign ha) $ cs++addHARule :: Rule ComponentSet+addHARule =+ describe "Rule for adding the HA component" .+ makeRule "component.ha" $ f+ where+ f :: ComponentSet -> Maybe ComponentSet+ f cs = do+ guard (derived cs `doesNotContain` AlternativeHypothesis)+ ha <- get AlternativeHypothesis (initials cs) + return $ append AlternativeHypothesis ha cs++----------------------------------------------------------+-- Rules for determining the properties of the data++chooseTTestRule :: Rule ComponentSet+chooseTTestRule =+ describe "Rule for choosing the type of test" .+ makeRule "component.test.t-test" $+ addTestChoice TTestOne++chooseTTestTwoRule :: Rule ComponentSet+chooseTTestTwoRule =+ describe "Rule for choosing the type of test" .+ makeRule "component.test.t-test-two" $+ addTestChoice TTestTwo++chooseTTestPairedRule :: Rule ComponentSet+chooseTTestPairedRule =+ describe "Rule for choosing the type of test" .+ makeRule "component.test.t-test-paired" $+ addTestChoice TTestPaired++-- | If the standard deviation of the population is known then the z-test is+-- /always/ chosen, thus a t-test is never chosen in that case. When this+-- information is unknown, then an estimation of the standard deviation of the+-- population has to be made based on the sample, in that case a t-test is+-- chosen.+--+-- When the sample size becomes large enough* it is also possible to choose a+-- z-test because the z-distribution looks like the t-distribution for large+-- sample sizes.+--+-- NOTE*: The threshold for what `large' means may vary. For now this threshold+-- is fixed at 100.+--+-- Consisely: If PopulationSdev known: ZTest+-- else if SampleSize big enough: ZTest or TTest+-- else TTest+chooseZTestRule :: Rule ComponentSet+chooseZTestRule =+ describe "Rule for choosing the type of test" .+ makeRule "component.test.z-test" $ f+ where+ f :: ComponentSet -> Maybe ComponentSet+ f cs =+ let largeThreshold = 100+ in do+ guard (derived cs `doesNotContain` TestChoice)+ n <- match naturalView $ fromMaybe 0 (getExpr SampleSize cs)+ if cs `contains` PopulationSdev || n >= largeThreshold + then return (append TestChoice (CChoice $ TestType ZTest) cs)+ else do+ -- This branch is similar to the previous rule body in revision 10549+ let tests = validTests cs+ guard (ZTest `elem` tests)+ return $ append TestChoice (CChoice $ TestType ZTest) cs++addTestFormulaRule :: Rule ComponentSet+addTestFormulaRule =+ describe "Rule for adding the test formula" .+ makeRule "component.test-formula" $ f+ where+ f :: ComponentSet -> [ComponentSet]+ f cs = do+ guard (derived cs `doesNotContain` TestFormula)+ rel <- inferTestFormula cs+ return $ append TestFormula (CRelation rel) cs++addTestValueRule :: Rule ComponentSet+addTestValueRule =+ describe "Rule for adding the test value (from the formula)" .+ makeRule "component.test-value" $ f+ where+ f :: ComponentSet -> [ComponentSet]+ f cs = do+ guard (derived cs `doesNotContain` TestValue)+ tv <- inferTestValue cs+ return $ append TestValue (CRelation tv) cs++inferTestValue :: MonadPlus m => ComponentSet -> m (Relation Expr)+inferTestValue cs =+ case getRelation TestValue (initials cs) of+ Just initialTestValue -> + return initialTestValue+ _ | inferTestChoice cs == Just ChiSquared -> + chiSquaredTestValue cs+ _ -> do+ new <- msum (map return $ applyAll addTestFormulaRule cs)+ let cs' = substitute new+ var <- leftHandSide <$> getRelation TestFormula new+ expr <- rightHandSide <$> getRelation TestFormula cs'+ val <- matchM doubleView expr+ return $ var .==. toExpr val++----------------------------------------------------------+-- Rules for performing a T-Test++addDfRule :: Rule ComponentSet+addDfRule =+ describe "Rule for adding the degrees of freedom component" .+ makeRule "component.df" $ f+ where+ f :: ComponentSet -> Maybe ComponentSet+ f cs = do+ guard (derived cs `doesNotContain` Df)+ case get Df (initials cs) of+ Just initialDf -> + return $ append Df initialDf cs+ _ | inferTestChoice cs == Just ChiSquared -> do+ df <- chiSquaredDf cs+ return $ append Df (CExpr (toExpr df)) cs+ _ -> do + test <- inferTestChoice cs+ df <- degreesOfFreedomFromTest test (chooseTypeOfTest cs)+ val <- matchM doubleView (getSubstitution cs |-> df)+ return $ append Df (CExpr (toExpr val)) cs++lookupTValueRule :: Rule ComponentSet+lookupTValueRule =+ describe "Rule for looking up a t-value" .+ makeRule "component.critical.t-value" $ + inferCriticalTWith $ \sided alpha df -> + [computeCriticalT sided alpha df]++-- shared function for computing t-value and r-value+inferCriticalGenericWith :: (TestType -> Bool) -> Expr -> (Sided -> Double -> Double -> [Double]) + -> ComponentSet -> [ComponentSet]+inferCriticalGenericWith forTestType var compute cs = do+ guard (cs `doesNotContain` Critical)+ guard (derived cs `contains` AlternativeHypothesis)+ guard (any forTestType (inferTestChoices cs))+ cs' <- matchM substitutedView cs+ alpha <- matchM doubleView <=< getExpr SignificanceLevel $ cs + let cs'' = case inferDf cs' of + Just df -> substitute (append Df (CExpr df) cs')+ Nothing -> cs'+ df <- matchM doubleView =<< getExpr Df cs''+ sided <- inferSidedness cs+ value <- compute sided alpha df+ return $ append Critical (CRelation $ var .==. fromDouble value) cs++inferCriticalTWith :: (Sided -> Double -> Double -> [Double]) -> ComponentSet -> [ComponentSet]+inferCriticalTWith = inferCriticalGenericWith isTTest (Var "tcrit")++computePValueTTest :: Rule ComponentSet+computePValueTTest =+ describe "Rule for computing the p-value for a t-test" .+ makeRule "component.p-value.t-test" $ f+ where+ f :: ComponentSet -> Maybe ComponentSet+ f cs = do+ guard (cs `doesNotContain` PValue)+ guard (derived cs `contains` AlternativeHypothesis)+ guard (maybe False isTTest (inferTestChoice cs))+ let cs' = substitute cs+ testStatistic <- match doubleView =<< fmap rightHandSide (getRelation TestValue cs')++ let cs'' = case inferDf cs' of + Just df -> substitute (append Df (CExpr df) cs')+ Nothing -> cs'+ df <- matchM doubleView =<< getExpr Df cs''+ sided <- inferSidedness cs+ let value = computePValueT sided testStatistic df+ return $ append PValue (CExpr $ fromDouble value) cs++----------------------------------------------------------+-- Rules for performing a Z-Test++lookupZValueRule :: Rule ComponentSet+lookupZValueRule =+ describe "Rule for looking up a z-value" .+ makeRule "component.critical.z-value" $ + inferCriticalZWith $ \sided alpha -> + [computeCriticalZ sided alpha]++inferCriticalZWith :: (Sided -> Double -> [Double]) -> ComponentSet -> [ComponentSet]+inferCriticalZWith compute cs = do+ guard (cs `doesNotContain` Critical)+ guard (derived cs `contains` AlternativeHypothesis)+ guard (maybe False (ZTest ==) (inferTestChoice cs))+ alpha <- matchM doubleView <=< getExpr SignificanceLevel $ cs+ sided <- inferSidedness cs+ value <- compute sided alpha+ return $ append Critical (CRelation $ Var "zcrit" .==. fromDouble value) cs++computePValueZTest :: Rule ComponentSet+computePValueZTest =+ describe "Rule for computing the p-value for a z-test" .+ makeRule "component.p-value.z-test" $ f+ where+ f :: ComponentSet -> Maybe ComponentSet+ f cs = do+ guard (cs `doesNotContain` PValue)+ guard (maybe False (ZTest ==) (inferTestChoice cs))+ guard (derived cs `contains` AlternativeHypothesis)+ let cs' = substitute cs+ -- TestStatistic was renamed to TestFormula, which is now a CRelation+ testStatistic <- match doubleView =<< fmap rightHandSide (getRelation TestValue cs')+ sided <- inferSidedness cs+ let value = computePValueZ sided testStatistic+ return $ append PValue (CExpr $ fromDouble value) cs++----------------------------------------------------------+-- Rules for performing R-Pearson++lookupRValueRule :: Rule ComponentSet+lookupRValueRule =+ describe "Rule for looking up a r-value" .+ makeRule "component.critical.r-value" $ + inferCriticalRWith $ \sided alpha df -> + [computeCriticalR sided alpha df]++chooseRPearsonRule :: Rule ComponentSet+chooseRPearsonRule =+ describe "Rule for choosing the type of test" .+ makeRule "component.test.r-pearson" $+ addTestChoice RPearson++addTestChoice :: TestType -> ComponentSet -> Maybe ComponentSet+addTestChoice testType cs = do+ guard (derived cs `doesNotContain` TestChoice)+ let tests = validTests cs+ guard (testType `elem` tests)+ return $ append TestChoice (CChoice $ TestType testType) cs++inferCriticalRWith :: (Sided -> Double -> Double -> [Double]) -> ComponentSet -> [ComponentSet]+inferCriticalRWith = inferCriticalGenericWith (== RPearson) (Var "rcrit")++----------------------------------------------------------+-- Rules for performing Anova++lookupFValueRule :: Rule ComponentSet+lookupFValueRule =+ describe "Rule for looking up a F-value" .+ makeRule "component.critical.f-value" $ inferCriticalFWith computeCriticalF++chooseAnovaRule :: Rule ComponentSet+chooseAnovaRule =+ describe "Rule for choosing the type of test" .+ makeRule "component.test.anova" $+ addTestChoice Anova++addDfBetweenWithin :: Rule ComponentSet+addDfBetweenWithin = describe "Add df between and within (for Anova)" $+ makeRule "component.df-anova" f+ where+ f :: ComponentSet -> Maybe ComponentSet+ f cs = do+ guard (derived cs `doesNotContain` DfBetween)+ guard (derived cs `doesNotContain` DfWithin)+ (between, within) <- inferDfBetweenWithin cs+ return $ append DfBetween (CExpr $ toExpr between) $ + append DfWithin (CExpr $ toExpr within) cs++inferDfBetweenWithin :: Monad m => ComponentSet -> m (Double, Double)+inferDfBetweenWithin cs = do+ n <- matchM doubleView <=< getExpr SampleSize $ cs+ let nrOfGroups = 2 -- always 2, for now+ dfBetween = nrOfGroups - 1+ dfWithin = n - nrOfGroups+ return (dfBetween, dfWithin)++inferCriticalFWith :: (Double -> Double -> Double -> [Double]) -> ComponentSet -> [ComponentSet]+inferCriticalFWith compute cs = do+ guard (cs `doesNotContain` Critical)+ guard (derived cs `contains` AlternativeHypothesis)+ guard (maybe False (Anova ==) (inferTestChoice cs))+ alpha <- matchM doubleView <=< getExpr SignificanceLevel $ cs+ (dfBetween, dfWithin) <- inferDfBetweenWithin cs+ value <- compute dfBetween dfWithin alpha+ return $ append Critical (CRelation $ Var "Fcrit" .==. fromDouble value) cs++----------------------------------------------------------+-- Rules for performing Chi-Squared++chooseChiSquaredRule :: Rule ComponentSet+chooseChiSquaredRule = describe "Rule for choosing the type of test" .+ makeRule "component.test.chi-squared" $+ addTestChoice ChiSquared++addHypothesesChiSquaredRule :: Rule ComponentSet+addHypothesesChiSquaredRule = describe "Add hypotheses (null and alternative) for chi-squared" $+ makeRule "component.hypotheses-chi-squared" f+ where+ f :: ComponentSet -> Maybe ComponentSet+ f cs = do+ guard (derived cs `doesNotContain` NullHypothesis)+ guard (derived cs `doesNotContain` AlternativeHypothesis)+ guard (maybe False (ChiSquared ==) (inferTestChoice cs))+ return $ append NullHypothesis (CExpr independent) + $ append AlternativeHypothesis (CExpr dependent) cs++addObservedTotals :: Rule ComponentSet+addObservedTotals = describe "Add totals (rows and columns) for observed frequencies" $ + makeRule "component.observed-totals" f+ where+ f :: ComponentSet -> Maybe ComponentSet+ f cs = do + guard (derived cs `doesNotContain` ObservedColumnTotals)+ guard (derived cs `doesNotContain` ObservedRowTotals)+ guard (derived cs `doesNotContain` ObservedTotal)+ table <- getTable ObservedFrequencies cs+ let (rowTotals, columnTotals, total) = computeTotals table+ return $ append ObservedRowTotals (CExpr $ toExpr rowTotals) + $ append ObservedColumnTotals (CExpr $ toExpr columnTotals) + $ append ObservedTotal (CExpr $ toExpr total) cs++type ChiSquaredTotals = ([Int], [Int], Int)++computeTotals :: [[Int]] -> ChiSquaredTotals+computeTotals table = + let rowTotals = map sum table+ columnTotals = map sum (transpose table)+ total = sum rowTotals+ in (rowTotals, columnTotals, total)++computeExpectedFrequencies :: ChiSquaredTotals -> [[Double]]+computeExpectedFrequencies (rowTotals, columnTotals, total) = table+ where+ xss = map (replicate (length columnTotals)) rowTotals+ yss = replicate (length rowTotals) columnTotals+ table = zipWith (zipWith f) xss yss+ f x y = fromIntegral (x*y) / fromIntegral total++addExpectedFrequencies :: Rule ComponentSet+addExpectedFrequencies = describe "Add expected frequencies" $ + makeRule "component.expected-frequencies" f+ where+ f :: ComponentSet -> Maybe ComponentSet+ f cs = do+ guard (derived cs `doesNotContain` ExpectedFrequencies)+ observed <- getTable ObservedFrequencies cs+ let totals = computeTotals observed+ expected = computeExpectedFrequencies totals+ return $ append ExpectedFrequencies (CExpr $ toExpr expected) cs+ +lookupChiValueRule :: Rule ComponentSet+lookupChiValueRule = describe "Rule for looking up a chi^2-value" .+ makeRule "component.critical.chi-value" $ + inferCriticalChiWith $ \sided alpha df -> + computeCriticalChi sided alpha df++inferCriticalChiWith :: (Sided -> Double -> Double -> [Double]) -> ComponentSet -> [ComponentSet]+inferCriticalChiWith = inferCriticalGenericWith (== ChiSquared) (Var "chicrit")++chiSquaredTestValue :: MonadPlus m => ComponentSet -> m (Relation Expr)+chiSquaredTestValue cs = do + observed <- getTable ObservedFrequencies cs+ expected <- getDoubleTable ExpectedFrequencies cs+ let table = zipWith (zipWith f) observed expected+ f o e = (fromIntegral o-e)^(2 :: Int) / e+ value = sum (map sum table)+ return (chiSquared .==. toExpr value)+ +chiSquaredDf :: MonadPlus m => ComponentSet -> m Int+chiSquaredDf cs = do+ observed <- getTable ObservedFrequencies cs+ guard (not $ null observed)+ let r = length observed+ c = length $ head observed+ return ((r-1)*(c-1))+ +getTable :: MonadPlus m => ComponentId -> ComponentSet -> m [[Int]]+getTable n cs = getExpr n cs >>= fromExpr++getDoubleTable :: MonadPlus m => ComponentId -> ComponentSet -> m [[Double]]+getDoubleTable n cs = getExpr n cs >>= fromExpr++----------------------------------------------------------+-- Rules for constructing the rejection region++-- lookup the Sidedness component; if not present, try to infer (using the rule)+inferSidedness :: MonadPlus m => ComponentSet -> m Sided+inferSidedness cs = + getSided Sidedness cs `mplus` do+ testType <- inferTestChoice cs+ guard (testType `elem` [Anova, ChiSquared])+ return RightSided+ `mplus` do+ cs' <- applyM determineSided cs+ getSided Sidedness cs'++inferVar :: MonadPlus m => ComponentSet -> m Expr+inferVar cs = do + testType <- inferTestChoice cs+ return (varForTestType testType)++varForTestType :: TestType -> Expr+varForTestType testType = + case testType of + ZTest -> Var "z"+ RPearson -> Var "r"+ Anova -> Var "F"+ ChiSquared -> chiSquared+ _ -> Var "t"++inferTestChoice :: MonadPlus m => ComponentSet -> m TestType+inferTestChoice cs =+ case validTests cs of+ hd:_ -> return hd+ _ -> fail "no valid test choice"++inferTestChoices :: ComponentSet -> [TestType]+inferTestChoices cs = + case getTestType TestChoice cs of+ Just test -> [test]+ Nothing -> validTests cs++inferTestFormula :: ComponentSet -> [Relation Expr]+inferTestFormula cs = + case getRelation TestFormula cs of+ Just formula -> [formula]+ Nothing -> do+ testType <- inferTestChoices cs+ let var = varForTestType testType+ let test = chooseTypeOfTest cs+ t <- testFormulaFromTest testType test+ return (var .==. t)++inferRejectionCritical :: MonadPlus m => ComponentSet -> m (Relation Expr)+inferRejectionCritical cs = + getRelation RejectionCritical cs `mplus` do+ cs' <- applyM addRejectionRule cs+ getRelation RejectionCritical cs'++inferConclusionCritical :: MonadPlus m => ComponentSet -> m Bool+inferConclusionCritical cs = + getConclusion ConclusionCritical cs `mplus` do+ cs' <- applyM criticalConclusionRule cs+ getConclusion ConclusionCritical cs'++inferConclusionPValue :: MonadPlus m => ComponentSet -> m (Relation Expr)+inferConclusionPValue cs = + getRelation ConclusionPValue cs `mplus` do+ cs' <- applyM addConclusionPValueRule cs+ getRelation ConclusionPValue cs'++inferDf :: MonadPlus m => ComponentSet -> m Expr+inferDf cs =+ getExpr Df cs `mplus` do+ cs' <- applyM addDfRule cs+ getExpr Df cs'++addRejectionRule :: Rule ComponentSet+addRejectionRule =+ describe "Rule for constructing the rejection critical component" .+ makeRule "component.rejection.critical" $ f+ where+ f :: ComponentSet -> [ComponentSet]+ f cs = do+ guard (cs `doesNotContain` RejectionCritical)+ guard (derived cs `contains` AlternativeHypothesis)+ let cs' = substitute cs+ sided <- inferSidedness cs'+ testType <- inferTestChoices cs+ let rel = case testType of + ZTest -> sidedRelation sided (Var "z") (Var "zcrit")+ RPearson -> sidedRelation sided (Var "r") (Var "rcrit")+ Anova -> sidedRelation sided (Var "F") (Var "Fcrit")+ ChiSquared -> sidedRelation sided chiSquared (Var "chicrit")+ _ -> sidedRelation sided (Var "t") (Var "tcrit")+ return . append RejectionCritical+ (CRelation rel) $ cs++addConclusionPValueRule :: Rule ComponentSet+addConclusionPValueRule =+ describe "Rule for constructing the conclusion p-value component" .+ makeRule "component.conclusion.p-value" $ f+ where+ f :: ComponentSet -> Maybe ComponentSet+ f cs = do+ guard (cs `doesNotContain` ConclusionPValue)+ pv <- matchM doubleView =<< getExpr PValue cs + alpha <- matchM doubleView =<< getExpr SignificanceLevel cs+ let relType = if pv > alpha then GreaterThan else LessThanOrEqualTo+ return . append ConclusionPValue+ (CRelation $ makeType relType (Var "p") (Var "alpha")) $ cs++----------------------------------------------------------+-- Rules for making a conclusion++criticalConclusionRule :: Rule ComponentSet+criticalConclusionRule =+ describe ("Rule for creating a conclusion based on a critical value " +++ "and test statistic") .+ makeRule "component.critical-conclusion" $ f+ where+ f :: ComponentSet -> Maybe ComponentSet+ f cs = do+ rej <- inferRejectionCritical cs+ -- to do: rejection critical is added to component set only to get the substituted relation+ let cs' = substitute (append RejectionCritical (CRelation rej) cs)+ rejection <- getRelation RejectionCritical cs'+ guard (cs `doesNotContain` ConclusionCritical)+ lhs <- match doubleView $ leftHandSide rejection+ rhs <- match doubleView $ rightHandSide rejection+ let result = eval (relationType rejection) lhs rhs+ return $ append ConclusionCritical (CChoice $ Conclusion result) cs++hypothesesConclusionCriticalRule :: Rule ComponentSet+hypothesesConclusionCriticalRule = + describe "derive the hypotheses conclusion from the critical conclusion" $+ makeRule "component.hypotheses-conclusion-critical" f + where+ f :: ComponentSet -> Maybe ComponentSet+ f cs = do+ guard (cs `doesNotContain` ConclusionHypotheses)+ concl <- inferConclusionCritical cs+ let rejhyp = if concl then RejectH0 else DontRejectH0+ return $ append ConclusionHypotheses (CChoice (RejectionHypotheses rejhyp)) cs++hypothesesConclusionPValueRule :: Rule ComponentSet+hypothesesConclusionPValueRule = + describe "derive the hypotheses conclusion from the p-value" $+ makeRule "component.hypotheses-conclusion-pvalue" f+ where+ f :: ComponentSet -> Maybe ComponentSet+ f cs = do+ guard (cs `doesNotContain` ConclusionHypotheses)+ rel <- inferConclusionPValue cs+ let rejhyp = if relationType rel == LessThanOrEqualTo then RejectH0 else DontRejectH0+ return $ append ConclusionHypotheses (CChoice (RejectionHypotheses rejhyp)) cs++addStandardErrorSigma :: Rule ComponentSet+addStandardErrorSigma = + describe "derive standard error from population standard deviation and sample size" $ + makeRule "component.standard-error-sigma" f + where+ f :: ComponentSet -> Maybe ComponentSet+ f cs = do+ guard (cs `doesNotContain` StandardError)+ n <- matchM doubleView =<< getExpr SampleSize cs+ psdev <- matchM doubleView =<< getRhsExpr PopulationSdev cs+ let se = Var "sigmaM" .==. toExpr (psdev / sqrt n)+ return $ append StandardError (CRelation se) cs+ +addStandardErrorSD :: Rule ComponentSet+addStandardErrorSD = + describe "derive standard error from sample standard deviation and sample size" $ + makeRule "component.standard-error-sd" f + where+ f :: ComponentSet -> Maybe ComponentSet+ f cs = do+ guard (cs `doesNotContain` StandardError)+ n <- matchM doubleView =<< getExpr SampleSize cs+ sdev <- matchM doubleView =<< getRhsExpr SampleSdev cs+ let se = Var "SEM" .==. toExpr (sdev / sqrt n)+ return $ append StandardError (CRelation se) cs
+ src/Domain/Hypothesis/Strategies.hs view
@@ -0,0 +1,74 @@+-----------------------------------------------------------------------------+-- Copyright 2020, Ideas project team. This file is distributed under the+-- terms of the Apache License 2.0. For more information, see the files+-- "LICENSE.txt" and "NOTICE.txt", which are included in the distribution.+-----------------------------------------------------------------------------+module Domain.Hypothesis.Strategies (hypothesisStrategy) where++import Domain.Hypothesis.Rules+import Domain.Statistics.Rules+import Domain.Statistics.ComponentSet+import Ideas.Common.Strategy.Combinators hiding (not, while)+import Ideas.Common.Library hiding (while)+import Prelude hiding (until, sequence, repeat)+++----------------------------------------------------------+-- Template for the hypothesis testing strategy++hypothesisStrategy :: LabeledStrategy ComponentSet+hypothesisStrategy = label "Hypothesis testing" $+ label "Preparation" (whileNotReady $ choice $+ [ addHypothesesRule, addH0FromHARule, addH0FromHAEqualSignRule, addHARule+ , addHypothesesChiSquaredRule+ , addAlphaRule, determineSided, chooseTTestRule+ , chooseTTestTwoRule, chooseTTestPairedRule, chooseZTestRule+ , chooseRPearsonRule, chooseAnovaRule, chooseChiSquaredRule+ ] ++ sampleStatistics)+ .*.+ check (\cs -> all (derived cs `contains`) [NullHypothesis, AlternativeHypothesis])+ .*.+ label "Computation" (whileNotReady $+ (check (\cs -> derived cs `doesNotContain` TestValue) .*. + (addTestFormulaRule .|. choice sampleStatistics))+ .|. (check allowCriticalRoute .*. choice+ [ addTestValueRule, addRejectionRule+ , lookupZValueRule, lookupTValueRule, lookupRValueRule, lookupFValueRule, lookupChiValueRule+ ])+ .|. (check allowPValueRoute .*. choice+ [ computePValueZTest, computePValueTTest+ ])+ )+ .*.+ check (\cs -> derived cs `contains` TestValue &&+ derived cs `contains` Critical || derived cs `contains` PValue)+ .*.+ label "Conclusion" (+ whileNotReady (criticalConclusionRule .|. addConclusionPValueRule)+ .*.+ (hypothesesConclusionCriticalRule .|. hypothesesConclusionPValueRule))+ where+ sampleStatistics = + [ addNRule, addAverageRule, addVarianceRule, addStandardDeviationRule+ , addStandardErrorRule, addStandardErrorSD, addDfRule, addDfBetweenWithin+ , addStandardErrorSigma, addObservedTotals, addExpectedFrequencies+ ]++-- customized while (combinator in library is greedy and uses repeat instead of many)+while :: IsStrategy f => (a -> Bool) -> f a -> Strategy a+while p s = many (check p .*. s)++whileNotReady :: Strategy ComponentSet -> Strategy ComponentSet+whileNotReady = while (`doesNotContain` ConclusionHypotheses)++allowCriticalRoute :: ComponentSet -> Bool+allowCriticalRoute cs = criticalRoute cs || not (pvalueRoute cs)++allowPValueRoute :: ComponentSet -> Bool+allowPValueRoute cs = pvalueRoute cs || not (criticalRoute cs)++criticalRoute :: ComponentSet -> Bool+criticalRoute cs = derived cs `contains` RejectionCritical || derived cs `contains` Critical++pvalueRoute :: ComponentSet -> Bool+pvalueRoute cs = derived cs `contains` PValue
+ src/Domain/Hypothesis/Tables.hs view
@@ -0,0 +1,88 @@+-----------------------------------------------------------------------------+-- Copyright 2020, Ideas project team. This file is distributed under the+-- terms of the Apache License 2.0. For more information, see the files+-- "LICENSE.txt" and "NOTICE.txt", which are included in the distribution.+-----------------------------------------------------------------------------+module Domain.Hypothesis.Tables+ ( zTableS, zTable, + tTableS, tTable,+ fTable, chiTable+ ) where++zTableS :: Bool -> Double -> Maybe Double+zTableS True alpha = zTable alpha+zTableS False alpha = zTable (alpha * 2)+ +zTable :: Double -> Maybe Double+zTable 0.05 = Just 1.645+zTable 0.025 = Just 1.960+zTable 0.01 = Just 2.326+zTable 0.005 = Just 2.576+zTable 0.0025 = Just 2.807+zTable 0.001 = Just 3.090+zTable _ = Nothing++tTableS :: Bool -> Double -> Int -> Maybe Double+tTableS True alpha df = tTable alpha df+tTableS False alpha df = tTable (alpha * 2) df++tTable :: Double -> Int -> Maybe Double+tTable alpha df | df > 60 = lookup alpha [(0.1, 1.282), (0.05, 1.645), (0.025, 1.960), (0.01, 2.326), (0.005, 2.576)]+ | otherwise = tTable' alpha df++tTable' :: Double -> Int -> Maybe Double+tTable' 0.1 df = lookup df [(1, 3.078), (2, 1.886), (3, 1.638), (4, 1.533), (5, 1.476), (6, 1.440), (7, 1.415), (8, 1.397), (9, 1.383), (10, 1.372),+ (11, 1.363), (12, 1.356), (13, 1.350), (14, 1.345), (15, 1.341), (16, 1.337), (17, 1.333), (18, 1.330), (19, 1.328), (20, 1.325),+ (21, 1.323), (22, 1.321), (23, 1.319), (24, 1.318), (25, 1.316), (26, 1.315), (27, 1.314), (28, 1.313), (29, 1.311), (30, 1.310),+ (40, 1.303), (60, 1.296)]+tTable' 0.05 df = lookup df [(1, 6.314), (2, 2.920), (3, 2.353), (4, 2.132), (5, 2.015), (6, 1.943), (7, 1.895), (8, 1.860), (9, 1.833), (10, 1.812),+ (11, 1.796), (12, 1.782), (13, 1.771), (14, 1.761), (15, 1.753), (16, 1.746), (17, 1.740), (18, 1.734), (19, 1.729), (20, 1.725),+ (21, 1.721), (22, 1.717), (23, 1.714), (24, 1.711), (25, 1.708), (26, 1.706), (27, 1.703), (28, 1.701), (29, 1.699), (30, 1.697),+ (40, 1.684), (60, 1.671)]+tTable' 0.025 df = lookup df [(1, 12.706), (2, 4.303), (3, 3.182), (4, 2.776), (5, 2.571), (6, 2.447), (7, 2.365), (8, 2.306), (9, 2.262), (10, 2.228),+ (11, 2.201), (12, 2.179), (13, 2.160), (14, 2.145), (15, 2.131), (16, 2.120), (17, 2.110), (18, 2.101), (19, 2.093), (20, 2.086),+ (21, 2.080), (22, 2.074), (23, 2.069), (24, 2.064), (25, 2.060), (26, 2.056), (27, 2.052), (28, 2.048), (29, 2.045), (30, 2.042),+ (40, 2.021), (60, 2.000)]+tTable' 0.01 df = lookup df [(1, 31.821), (2, 6.965), (3, 4.541), (4, 3.747), (5, 3.365), (6, 3.143), (7, 2.998), (8, 2.896), (9, 2.821), (10, 2.764),+ (11, 2.718), (12, 2.681), (13, 2.650), (14, 2.624), (15, 2.602), (16, 2.583), (17, 2.567), (18, 2.552), (19, 2.539), (20, 2.528),+ (21, 2.518), (22, 2.508), (23, 2.500), (24, 2.492), (25, 2.485), (26, 2.479), (27, 2.473), (28, 2.467), (29, 2.462), (30, 2.457),+ (40, 2.423), (60, 2.390)]+tTable' 0.005 df = lookup df [(1, 63.657), (2, 9.925), (3, 5.841), (4, 4.604), (5, 4.032), (6, 3.707), (7, 3.499), (8, 3.355), (9, 3.250), (10, 3.169),+ (11, 3.106), (12, 3.055), (13, 3.012), (14, 2.997), (15, 2.947), (16, 2.921), (17, 2.898), (18, 2.878), (19, 2.861), (20, 2.845),+ (21, 2.831), (22, 2.819), (23, 2.807), (24, 2.797), (25, 2.787), (26, 2.779), (27, 2.771), (28, 2.763), (29, 2.756), (30, 2.750),+ (40, 2.704), (60, 2.660)]+tTable' _ _ = Nothing++-- for Anova+fTable :: Double -> Double -> Double -> Maybe Double+fTable dfBetween dfWithin alpha+ | alpha == 0.05 = lookup (dfBetween, dfWithin) + [((1, 38), 4.10), ((1, 40), 4.08), ((2, 38), 3.25), ((2, 40), 3.23) ]+ | alpha == 0.01 = lookup (dfBetween, dfWithin) + [((1, 38), 7.35), ((1, 40), 7.31), ((2, 38), 5.21), ((2, 40), 5.18) ]+ | otherwise = Nothing+ +-- for Chi squared+chiTable :: Double -> Int -> Maybe Double+chiTable alpha df+ | alpha == 0.10 && df > 0 && df <= 20 = Just $+ [ 2.71, 4.61, 6.25, 7.78, 9.24, 10.64, 12.02, 13.36, 14.68, 15.99 + , 17.28, 18.55, 19.81, 21.06, 22.31, 23.54, 24.77, 25.99, 27.20, 28.41+ ] !! (df - 1)+ | alpha == 0.05 && df > 0 && df <= 20 = Just $+ [ 3.84, 5.99, 7.81, 9.49, 11.07, 12.59, 14.07, 15.51, 16.92, 18.31+ , 19.68, 21.03, 22.36, 23.68, 25.00, 26.30, 27.59, 28.87, 30.14, 31.41+ ] !! (df - 1)+ | alpha == 0.025 && df > 0 && df <= 20 = Just $+ [ 5.02, 7.38, 9.35, 11.14, 12.83, 14.45, 16.01, 17.53, 19.02, 20.48+ , 21.92, 23.34, 24.74, 26.12, 27.49, 28.85, 30.19, 31.53, 32.85, 34.17+ ] !! (df - 1)+ | alpha == 0.01 && df > 0 && df <= 20 = Just $+ [ 6.63, 9.21, 11.34, 13.28, 15.09, 16.81, 18.48, 20.09, 21.67, 23.21+ , 24.72, 26.22, 27.69, 29.14, 30.58, 32.00, 33.41, 34.81, 36.19, 37.57+ ] !! (df - 1)+ | alpha == 0.005 && df > 0 && df <= 20 = Just $+ [ 7.88, 10.60, 12.84, 14.86, 16.75, 18.55, 20.28, 21.96, 23.59, 25.19+ , 26.76, 28.30, 29.82, 31.32, 32.80, 34.27, 35.72, 37.16, 38.58, 40.00+ ] !! (df - 1)+ | otherwise = Nothing
+ src/Domain/Statistics/Component.hs view
@@ -0,0 +1,328 @@+-----------------------------------------------------------------------------+-- Copyright 2020, Ideas project team. This file is distributed under the+-- terms of the Apache License 2.0. For more information, see the files+-- "LICENSE.txt" and "NOTICE.txt", which are included in the distribution.+-----------------------------------------------------------------------------+{-# LANGUAGE FlexibleInstances #-}++module Domain.Statistics.Component+ ( Choice (..), Component (..), ComponentId (..)+ , Sided (..), TestType (..), RejectionHypotheses(..)+ , choices, getCId, fromCId, sidedRelation+ , isData, isExpr, isRelation+ , isTestType, isSided, isRejectionHypotheses, isConclusion+ ) where++import Control.Applicative ((<|>))+import Control.Monad+import Data.Char+import Data.Function+import Data.List+import Domain.Statistics.Symbols+import Domain.Math.Data.Relation+import Domain.Math.Expr+import Domain.Statistics.Data+import Ideas.Common.Rewriting hiding (trueSymbol,falseSymbol)+import Prelude hiding (lookup)++data Component+ = CData Data -- Data component type (sample data)+ | CExpr Expr -- Expression component type+ | CRelation (Relation Expr)+ | CChoice Choice -- Expressing various choices+ deriving Eq++instance Show Component where+ show (CRelation r) = show r+ show (CExpr e) = show e+ show (CChoice c) = show c+ show (CData d) = show d++instance IsTerm Component where+ toTerm (CData x) = unary cDataSymbol (toTerm x)+ toTerm (CChoice x) = toTerm x+ toTerm (CRelation x) = toTerm x+ toTerm (CExpr x) = toTerm x++ -- Check to see if the whole x can be parsed as an expression+ fromTerm (TCon s [val])+ | s == cDataSymbol = CData <$> fromTerm val+ -- First try to convert it to a relation because relation is more specific+ -- than an Expression (IMPORTANT)+ fromTerm t = CChoice <$> fromTerm t <|>+ CRelation <$> fromTerm t <|>+ CExpr <$> fromTerm t++data Choice+ = TestType TestType+ | Sided Sided+ | RejectionHypotheses RejectionHypotheses+ | Conclusion Bool+ deriving Eq++instance Show Choice where+ show (TestType tt) = show tt+ show (Sided sd) = show sd+ show (RejectionHypotheses rh) = show rh+ show (Conclusion c) = map toLower (show c)++instance IsTerm Choice where+ toTerm (TestType x) = toTerm x+ toTerm (Sided x) = toTerm x+ toTerm (RejectionHypotheses rh) = toTerm rh+ toTerm (Conclusion x) = toTermBool x++ fromTerm t =+ Sided <$> fromTerm t <|>+ TestType <$> fromTerm t <|>+ RejectionHypotheses <$> fromTerm t <|>+ Conclusion <$> fromTermBool t++choices :: [Choice]+choices = [ TestType tt | tt <- [minBound .. maxBound] ] +++ [ Sided sd | sd <- [minBound .. maxBound] ] +++ [ RejectionHypotheses rj | rj <- [minBound .. maxBound] ] +++ [ Conclusion c | c <- [minBound .. maxBound] ]++data TestType+ = TTestOne+ | TTestTwo+ | TTestPaired+ | ZTest+ | RPearson+ | Anova+ | ChiSquared+ deriving (Eq, Enum, Bounded)++instance Show TestType where+ -- to do: use same string as corresponding symbols (without dashes)+ show TTestOne = "t-test-one-sample"+ show TTestTwo = "t-test-two-sample"+ show TTestPaired = "t-test-paired"+ show ZTest = "z-test"+ show RPearson = "r-pearson"+ show Anova = "anova"+ show ChiSquared = "chi-squared"++data Sided+ = TwoSided+ | LeftSided+ | RightSided+ deriving (Eq, Enum, Bounded)++instance Show Sided where+ show TwoSided = "two"+ show LeftSided = "left"+ show RightSided = "right"++sidedRelation :: Sided -> Expr -> Expr -> Relation Expr+sidedRelation LeftSided lhs rhs = makeType LessThan lhs rhs+sidedRelation RightSided lhs rhs = makeType GreaterThan lhs rhs+sidedRelation TwoSided lhs rhs = makeType GreaterThan (abs lhs) rhs++data RejectionHypotheses+ = RejectH0+ | AcceptH0+ | DontRejectH0+ | RejectH1+ | AcceptH1+ | DontRejectH1+ deriving (Eq, Enum, Bounded)+ +instance Show RejectionHypotheses where+ show RejectH0 = "rejecth0" + show AcceptH0 = "accepth0" + show DontRejectH0 = "dontrejecth0" + show RejectH1 = "rejecth1" + show AcceptH1 = "accepth1" + show DontRejectH1 = "dontrejecth1"++data ComponentId+ = DataSet+ | SampleSize -- n+ | SampleMean -- M+ | SampleVariance -- s^2+ | SampleSdev -- s+ | ObservedFrequencies -- fo (two-dimensional matrix), for chi squared+ | ObservedRowTotals+ | ObservedColumnTotals+ | ObservedTotal+ | ExpectedFrequencies -- fe (two-dimensional matrix), for chi squared+ | PopulationMean -- mu+ | PopulationSdev -- sigma+ | StandardError -- sigmaM or SEM+ | SignificanceLevel -- alpha+ | NullHypothesis+ | AlternativeHypothesis+ | PooledVariance -- S_p+ | Df -- df+ | DfBetween -- for Anova+ | DfWithin -- for Anova+ | Groups+ | Sidedness+ | TestChoice+ | Critical+ | PValue+ | RejectionCritical+ | ConclusionPValue+ | ConclusionCritical+ | ConclusionHypotheses+ | TestFormula+ | TestValue+ | Proportion -- P+ | Correlation -- r+ | One ComponentId -- Belonging to the first sample (in case of multiple samples)+ | Two ComponentId -- Belonging to the second sample+ | Other String+ deriving (Eq, Show)++instance Ord ComponentId where+ compare = compare `on` getCId++instance IsTerm ComponentId where+ toTerm = toTerm . getCId+ fromTerm t = fromCId <$> fromTerm t++getCId :: ComponentId -> String+getCId (Other nm) = nm+getCId (One nm) = getCId nm ++ "1"+getCId (Two nm) = getCId nm ++ "2"+getCId cId = case lookup cId cIdTable of+ Just nm -> nm+ _ -> error "Missing name for ComponentId"+fromCId :: String -> ComponentId+fromCId s+ | suffix == "1" = One (fromCId $ init s)+ | suffix == "2" = Two (fromCId $ init s)+ | otherwise =+ case [ x | (x, y) <- cIdTable, y == s ] of+ [x] -> x+ _ -> Other s+ where+ suffix = take 1 (reverse s)++cIdTable :: [(ComponentId, String)]+cIdTable =+ [ (DataSet, "data")+ , (SampleSize, "n")+ , (SampleMean, "samplemean") -- was: "mean"+ , (SampleVariance, "samplevariance") -- was: "variance"+ , (SampleSdev, "samplesd") -- was: "sd"+ , (ObservedFrequencies, "fo")+ , (ObservedRowTotals, "forowtotal")+ , (ObservedColumnTotals, "focolumntotal")+ , (ObservedTotal, "fototal")+ , (ExpectedFrequencies, "fe")+ , (PopulationMean, "populationmean") -- was: "pmean"+ , (PopulationSdev, "populationsd") -- was "psd"+ , (StandardError, "standarderror") -- was: "se"+ , (SignificanceLevel, "alpha")+ , (NullHypothesis, "h0")+ , (AlternativeHypothesis, "ha")+ , (PooledVariance, "pooledvariance") -- was "sp2"+ , (Df, "df")+ , (DfWithin, "dfwithin")+ , (DfBetween, "dfbetween")+ , (Groups, "groups")+ , (Sidedness, "sided")+ , (TestChoice, "test")+ , (Critical, "critical")+ , (PValue, "pvalue")+ , (RejectionCritical, "rejectioncritical") -- was: "rejection"+ , (ConclusionPValue, "conclusionpvalue")+ , (ConclusionCritical, "conclusioncritical") -- was: "conclusion"+ , (ConclusionHypotheses, "conclusionhypotheses")+ , (TestFormula, "testformula")+ , (TestValue, "testvalue")+ , (Proportion, "p")+ , (Correlation, "r")+ ]++----------------------------------------------------------+-- Component utilities++isData :: Monad m => Component -> m Data+isData (CData d) = return d+isData _ = failMsg "data"++isExpr :: Monad m => Component -> m Expr+isExpr (CExpr e) = return e+isExpr _ = failMsg "expr"++isRelation :: Monad m => Component -> m (Relation Expr)+isRelation (CRelation r) = return r+isRelation _ = failMsg "relation"++isTestType :: Monad m => Component -> m TestType+isTestType (CChoice (TestType tt)) = return tt+isTestType _ = failMsg "testtype"++isSided :: Monad m => Component -> m Sided+isSided (CChoice (Sided s)) = return s+isSided _ = failMsg "sided"++isRejectionHypotheses :: Monad m => Component -> m RejectionHypotheses+isRejectionHypotheses (CChoice (RejectionHypotheses h)) = return h+isRejectionHypotheses _ = failMsg "rejectionhypotheses"++isConclusion :: Monad m => Component -> m Bool+isConclusion (CChoice (Conclusion b)) = return b+isConclusion _ = failMsg "conclusion"++failMsg :: Monad m => String -> m a+failMsg s = fail ("component is not of type " ++ s)++----------------------------------------------------------+-- Component utilities++-- | Types that have only atomic constructors are easy and follow a pattern for+-- their implementation+instance IsTerm TestType where+ toTerm TTestOne = symbol tTestOneSymbol+ toTerm TTestTwo = symbol tTestTwoSymbol+ toTerm TTestPaired = symbol tTestPairedSymbol+ toTerm ZTest = symbol zTestSymbol+ toTerm RPearson = symbol rPearsonSymbol+ toTerm Anova = symbol anovaSymbol+ toTerm ChiSquared = symbol chiSquaredSymbol++ fromTerm (TCon s [])+ | s == tTestOneSymbol = return TTestOne+ | s == tTestTwoSymbol = return TTestTwo+ | s == tTestPairedSymbol = return TTestPaired+ | s == zTestSymbol = return ZTest+ | s == rPearsonSymbol = return RPearson+ | s == anovaSymbol = return Anova+ | s == chiSquaredSymbol = return ChiSquared+ fromTerm (TVar "rpearson") = return RPearson -- temporary fix!!!+ fromTerm (TVar "anova") = return Anova -- temporary fix!!!+ fromTerm (TVar "chisquared") = return ChiSquared -- temporary fix!!!+ fromTerm t = fail $ "Invalid term: " ++ show t ++ " not TestType"++instance IsTerm Sided where+ toTerm TwoSided = symbol twoSidedSymbol+ toTerm LeftSided = symbol leftSidedSymbol+ toTerm RightSided = symbol rightSidedSymbol++ fromTerm (TCon s [])+ | s == twoSidedSymbol = return TwoSided+ | s == leftSidedSymbol = return LeftSided+ | s == rightSidedSymbol = return RightSided+ fromTerm t = fail $ "Invalid term: " ++ show t ++ " not Sided"+ +instance IsTerm RejectionHypotheses where+ toTerm RejectH0 = symbol rejectH0Symbol+ toTerm AcceptH0 = symbol acceptH0Symbol+ toTerm DontRejectH0 = symbol dontRejectH0Symbol+ toTerm RejectH1 = symbol rejectH1Symbol+ toTerm AcceptH1 = symbol acceptH1Symbol+ toTerm DontRejectH1 = symbol dontRejectH1Symbol+ + fromTerm (TCon s [])+ | s == rejectH0Symbol = return RejectH0+ | s == acceptH0Symbol = return AcceptH0+ | s == dontRejectH0Symbol = return DontRejectH0+ | s == rejectH1Symbol = return RejectH1+ | s == acceptH1Symbol = return AcceptH1+ | s == dontRejectH1Symbol = return DontRejectH1+ fromTerm t = fail $ "Invalid term: " ++ show t ++ " not RejectionHypotheses"
+ src/Domain/Statistics/ComponentSet.hs view
@@ -0,0 +1,163 @@+-----------------------------------------------------------------------------+-- Copyright 2020, Ideas project team. This file is distributed under the+-- terms of the Apache License 2.0. For more information, see the files+-- "LICENSE.txt" and "NOTICE.txt", which are included in the distribution.+-----------------------------------------------------------------------------+module Domain.Statistics.ComponentSet+ ( -- * Types+ ComponentSet, module Domain.Statistics.Component+ -- * Constructors+ , initialSet, derivedSet+ -- * Combining and transforming+ , append, delete, mapComponent, initials, derived, toList+ -- * Getters+ , get, getData, getExpr, getRelation, getRhsExpr+ , getTestType, getSided, getRejectionHypotheses, getConclusion+ -- * Membership+ , contains, doesNotContain+ , isInitial, isDerived+ ) where++import Control.Monad+import Data.List hiding (delete)+import Data.Maybe+import Domain.Statistics.Symbols+import Domain.Math.Data.Relation+import Domain.Math.Expr+import Domain.Statistics.Data+import Domain.Statistics.Component+import Ideas.Common.Rewriting hiding (trueSymbol,falseSymbol)+import Ideas.Common.Classes++------------------------------------------------------------------------------+-- Types++data ComponentSet = CS + { initialsList :: [(ComponentId, Component)]+ , derivedList :: [(ComponentId, Component)]+ }++instance Show ComponentSet where+ show cs = unlines (map (f True) (initialsList cs) ++ map (f False) (derivedList cs))+ where+ f b (k, v) = bracketsIf b (getCId k) ++ ": " ++ show v ++ ";"++ bracketsIf True s = "[" ++ s ++ "]"+ bracketsIf False s = s++instance Eq ComponentSet where+ CS xs1 xs2 == CS ys1 ys2 = f xs1 == f ys1 && f xs2 == f ys2+ where+ f = sortBy (\(x, _) (y, _) -> x `compare` y)++instance Semigroup ComponentSet where+ CS xs1 xs2 <> CS ys1 ys2 = CS (xs1 ++ ys1) (xs2 ++ ys2)++instance Monoid ComponentSet where+ mempty = CS [] []+ mappend = (<>)++instance IsTerm ComponentSet where+ toTerm cs = TList $ + [ ternary componentSymbol (toTerm k) (symbol initialSymbol) (toTerm c) + | (k, c) <- initialsList cs+ ] +++ [ ternary componentSymbol (toTerm k) (symbol derivedSymbol) (toTerm c) + | (k, c) <- derivedList cs+ ]++ fromTerm (TList xs) =+ let f (TCon s [t1, TCon sv [], t3]) | s == componentSymbol = do+ k <- fromTerm t1+ c <- fromTerm t3+ if sv == initialSymbol + then return $ CS [(k, c)] []+ else if sv == derivedSymbol + then return $ CS [] [(k, c)]+ else fail "expected initial/derived"+ f _ = fail "Term is not a component"+ in mconcat <$> mapM f xs+ fromTerm _ = fail "Term is not a component (not a list)"++------------------------------------------------------------------------------+-- Constructors++initialSet :: [(ComponentId, Component)] -> ComponentSet+initialSet xs = CS xs []++derivedSet :: [(ComponentId, Component)] -> ComponentSet+derivedSet = CS []++------------------------------------------------------------------------------+-- Combining and transforming++append :: ComponentId -> Component -> ComponentSet -> ComponentSet+append k v cs = cs <> derivedSet [(k, v)]++delete :: ComponentId -> ComponentSet -> ComponentSet+delete k (CS xs ys) = CS (f xs) (f ys)+ where+ f = filter ((/= k) . fst)++mapComponent :: (Component -> Component) -> ComponentSet -> ComponentSet+mapComponent f (CS xs ys) = CS (change xs) (change ys)+ where+ change = map (mapSecond f)+ +initials :: ComponentSet -> ComponentSet+initials cs = mempty { initialsList = initialsList cs }++derived :: ComponentSet -> ComponentSet +derived cs = mempty { derivedList = derivedList cs }++toList :: ComponentSet -> [(ComponentId, Component)]+toList cs = initialsList cs ++ derivedList cs++------------------------------------------------------------------------------+-- Getters++get :: Monad m => ComponentId -> ComponentSet -> m Component+get n cs = + -- first consider derived components+ case lookup n (derivedList cs ++ initialsList cs) of+ Just c -> return c+ Nothing -> fail $ show n ++ " missing"++getData :: Monad m => ComponentId -> ComponentSet -> m Data+getData n = get n >=> isData++getExpr :: Monad m => ComponentId -> ComponentSet -> m Expr+getExpr n = get n >=> isExpr++getRelation :: Monad m => ComponentId -> ComponentSet -> m (Relation Expr)+getRelation n = get n >=> isRelation++getRhsExpr :: MonadPlus m => ComponentId -> ComponentSet -> m Expr+getRhsExpr n cs = fmap rightHandSide (getRelation n cs) `mplus` getExpr n cs++getTestType :: Monad m => ComponentId -> ComponentSet -> m TestType+getTestType n = get n >=> isTestType++getSided :: Monad m => ComponentId -> ComponentSet -> m Sided+getSided n = get n >=> isSided++getConclusion :: Monad m => ComponentId -> ComponentSet -> m Bool+getConclusion n = get n >=> isConclusion++getRejectionHypotheses :: Monad m => ComponentId -> ComponentSet -> m RejectionHypotheses+getRejectionHypotheses n = get n >=> isRejectionHypotheses++------------------------------------------------------------------------------+-- Membership++contains :: ComponentSet -> ComponentId -> Bool+contains cs x = isJust (get x cs)++doesNotContain :: ComponentSet -> ComponentId -> Bool+doesNotContain cs = not . contains cs++isInitial :: ComponentId -> ComponentSet -> Bool+isInitial n cs = initials cs `contains` n++isDerived :: ComponentId -> ComponentSet -> Bool+isDerived n cs = derived cs `contains` n
+ src/Domain/Statistics/Data.hs view
@@ -0,0 +1,86 @@+-----------------------------------------------------------------------------+-- Copyright 2020, Ideas project team. This file is distributed under the+-- terms of the Apache License 2.0. For more information, see the files+-- "LICENSE.txt" and "NOTICE.txt", which are included in the distribution.+-----------------------------------------------------------------------------+module Domain.Statistics.Data+ ( Data (..)+ , isSingleSample+ , sampleSize+ , sampleSum+ , sampleMean+ , sampleVariance+ ) where++import Ideas.Common.Rewriting+import Domain.Statistics.Symbols++-- | Data component+data Data+ = SingleSample [Double] -- Single sample+ | TwoSample [Double] [Double] Bool -- Two samples, either paired or not+ deriving (Eq, Show, Read)++instance IsTerm Data where+ toTerm (SingleSample xs) = unary singleSampleSymbol (toTerm xs)+ toTerm (TwoSample xs ys b) =+ ternary twoSampleSymbol (toTerm xs) (toTerm ys) (toTermBool b)++ fromTerm (TCon s [x])+ | s == singleSampleSymbol = + SingleSample <$> fromTerm x+ fromTerm (TCon s [xs, ys, b]) + | s == twoSampleSymbol = + TwoSample <$> fromTerm xs <*> fromTerm ys <*> fromTermBool b+ fromTerm t = fail $ "Invalid term: " ++ show t ++ " not Data"+ +-- | Utility functions+isSingleSample :: Data -> Bool+isSingleSample (SingleSample _) = True+isSingleSample _ = False++-- |+sampleSize :: Data -> Int+sampleSize (SingleSample xs) = length xs+sampleSize (TwoSample xs _ True) = length xs -- Both samples should have the same sample size+sampleSize _ = error "two sample sizes for unpaired samples"++-- |+sampleSum :: Data -> [Double]+sampleSum (SingleSample xs) = [sum xs]+sampleSum (TwoSample xs ys _) = [sum xs, sum ys]++-- |+sampleMean :: Data -> [Double]+sampleMean (SingleSample xs) = [ sum xs / fromIntegral (length xs) ]+sampleMean (TwoSample xs ys _) = [ sum xs / fromIntegral (length xs)+ , sum ys / fromIntegral (length ys)+ ]++-- |+sampleVariance :: Bool -> Data -> [Double]+sampleVariance True = sampleVarianceT+sampleVariance False = sampleVarianceZ++-- |+sampleVarianceZ :: Data -> [Double]+sampleVarianceZ (SingleSample xs) =+ let mean = head $ sampleMean (SingleSample xs)+ in [foldr (\x r -> r + (x - mean) ** 2) 0 xs / fromIntegral (length xs)]+sampleVarianceZ (TwoSample xs ys _) =+ let (mean1:mean2:_) = sampleMean (TwoSample xs ys False)+ in [ foldr (\x r -> r + (x - mean1) ** 2) 0 xs / fromIntegral (length xs)+ , foldr (\x r -> r + (x - mean2) ** 2) 0 ys / fromIntegral (length ys)+ ]++-- |+sampleVarianceT :: Data -> [Double]+sampleVarianceT (SingleSample xs) =+ let mean = head $ sampleMean (SingleSample xs)+ in [foldr (\x r -> r + (x - mean) ** 2) 0 xs / fromIntegral (length xs - 1)]+sampleVarianceT (TwoSample xs ys _) =+ let (mean1:mean2:_) = sampleMean (TwoSample xs ys False)+ in [ foldr (\x r -> r + (x - mean1) ** 2) 0 xs / fromIntegral (length xs - 1)+ , foldr (\x r -> r + (x - mean2) ** 2) 0 ys / fromIntegral (length ys - 1)+ ]+
+ src/Domain/Statistics/Parser.hs view
@@ -0,0 +1,71 @@+-----------------------------------------------------------------------------+-- Copyright 2020, Ideas project team. This file is distributed under the+-- terms of the Apache License 2.0. For more information, see the files+-- "LICENSE.txt" and "NOTICE.txt", which are included in the distribution.+-----------------------------------------------------------------------------+module Domain.Statistics.Parser (parseComponentSet) where++import Domain.Math.Expr+import Domain.Statistics.ComponentSet+import Ideas.Utils.Parsing hiding (char)+import qualified Ideas.Utils.Parsing as P++parseComponentSet :: String -> Either String ComponentSet+parseComponentSet = parseSimple pComponentSet ++pComponentSet :: Parser ComponentSet+pComponentSet = mconcat <$> many pComponent++pComponent :: Parser ComponentSet+pComponent = do+ (k, b) <- pKey+ char ':'+ v <- pValue+ char ';'+ return $ if b then initialSet [(k, v)] else derivedSet [(k, v)]++pKey :: Parser (ComponentId, Bool)+pKey = (\n -> (n, True)) <$> brackets pComponentId + <|> (\n -> (n, False)) <$> pComponentId++pComponentId :: Parser ComponentId+pComponentId = fromCId <$> parseName++pValue :: Parser Component+pValue = CChoice <$> try parseChoice+ <|> CRelation <$> try pRelExpr+ <|> CExpr <$> try pExpr+ +brackets :: Parser a -> Parser a+brackets p = do+ char '['+ a <- p+ char ']'+ return a++parseChoice :: Parser Choice+parseChoice = choice [ x <$ try (string (show x)) | x <- choices ]++parseName :: Parser String+parseName =+ let+ firstChar = oneOf $ ['a'..'z'] ++ ['A'..'Z']+ nonFirstChar = oneOf $ ['a'..'z'] ++ ['A'..'Z'] ++ ['0'..'9'] ++ "-"+ in do+ fc <- firstChar+ rest <- many nonFirstChar+ return (fc:rest)++----------------------------------------------------------------------------++char :: Char -> Parser ()+char c = lexeme (P.char c) >> return ()++whitespace :: Parser ()+whitespace = many (oneOf " \n\t") >> return ()++lexeme :: Parser a -> Parser a+lexeme p = do+ x <- p+ whitespace+ return x
+ src/Domain/Statistics/Rules.hs view
@@ -0,0 +1,94 @@+-----------------------------------------------------------------------------+-- Copyright 2020, Ideas project team. This file is distributed under the+-- terms of the Apache License 2.0. For more information, see the files+-- "LICENSE.txt" and "NOTICE.txt", which are included in the distribution.+-----------------------------------------------------------------------------+module Domain.Statistics.Rules+ ( addNRule+ , addAverageRule+ , addVarianceRule+ , addStandardDeviationRule+ , addStandardErrorRule+ ) where++import Control.Monad+import Domain.Math.Expr+import Domain.Statistics.ComponentSet+import Ideas.Common.Library+import Domain.Statistics.Data++----------------------------------------------------------+-- Rules on the data component++addNRule :: Rule ComponentSet+addNRule = describe "Rule for computing the n of the dataset" $ makeRule "component.n" f+ where+ f :: ComponentSet -> Maybe ComponentSet+ f cs = do guard (not $ contains cs SampleSize)+ sample <- getData DataSet cs+ let nCs = nComponents sample+ return $ foldr (\(nm, c) cs' -> append nm c cs') cs nCs+ nComponents (SingleSample xs) = [(SampleSize, CExpr $ Nat $ toInteger (length xs))]+ nComponents (TwoSample xs ys _) = [(One SampleSize, CExpr $ Nat $ toInteger (length xs)),+ (Two SampleSize, CExpr $ Nat $ toInteger (length ys))]++addAverageRule :: Rule ComponentSet+addAverageRule = describe "Rule for adding the formula for computing the mean" $ makeRule "component.mean" f+ where+ f :: ComponentSet -> Maybe ComponentSet+ f cs = do guard (not $ contains cs SampleMean)+ guard (not $ contains cs (One SampleMean))+ guard (not $ contains cs (Two SampleMean))+ sample <- getData DataSet cs+ -- n <- getExpr SampleMean cs+ case sampleMean sample of+ [mean] -> return $ append SampleMean (CExpr $ fromDouble mean) cs+ [mean0, mean1] -> return $ append (Two SampleMean) (CExpr $ fromDouble mean1) (append (One SampleMean) (CExpr $ fromDouble mean0) cs)+ _ -> Nothing+++addVarianceRule :: Rule ComponentSet+addVarianceRule = describe "Rule for adding the formula for computing the variance" $ makeRule "component.variance" f+ where+ f :: ComponentSet -> Maybe ComponentSet+ f cs = do guard (not $ contains cs SampleVariance)+ guard (not $ contains cs (One SampleVariance))+ guard (not $ contains cs (Two SampleVariance))+ sample <- getData DataSet cs+ -- n <- getExpr SampleSize cs+ -- mean <- getExpr SampleMean cs+ test <- getTestType TestChoice cs+ -- s = sum((X_i - avg(X))^2) / n+ case sampleVariance (test /= ZTest) sample of+ [var] -> return $ append SampleVariance (CExpr $ fromDouble var) cs+ [var0, var1] -> return $ append (Two SampleVariance) (CExpr $ fromDouble var1) (append (One SampleVariance) (CExpr $ fromDouble var0) cs)+ _ -> Nothing++addStandardDeviationRule :: Rule ComponentSet+addStandardDeviationRule = describe "Rule for adding the formula for computing the deviation" $ makeRule "component.standard-deviation" f+ where+ f :: ComponentSet -> Maybe ComponentSet+ f cs = do guard (not $ contains cs SampleSdev)+ guard (not $ contains cs (One SampleSdev))+ guard (not $ contains cs (Two SampleSdev))+ sample <- getData DataSet cs+ -- n <- getExpr SampleSize cs+ -- mean <- getExpr SampleMean cs+ test <- getTestType TestChoice cs+ -- s = sqrt(sum((X_i - avg(X))^2) / n)+ case map sqrt $ sampleVariance (test /= ZTest) sample of+ [sd] -> return $ append SampleSdev (CExpr $ fromDouble sd) cs+ [sd0, sd1] -> return $ append (Two SampleSdev) (CExpr $ fromDouble sd1) (append (One SampleSdev) (CExpr $ fromDouble sd0) cs)+ _ -> Nothing++addStandardErrorRule :: Rule ComponentSet+addStandardErrorRule = describe "Rule for adding the formula for computing the standard error" $ makeRule "component.standard-error" f+ where+ f :: ComponentSet -> Maybe ComponentSet+ f cs = do guard (not $ contains cs StandardError)+ guard (not $ contains cs (One StandardError))+ guard (not $ contains cs (Two StandardError))+ sd <- getExpr SampleSdev cs+ -- n <- getExpr SampleSize cs+ -- SE = s / sqrt(n)+ return $ append StandardError (CExpr $ sd / sqrt(Var "n")) cs
+ src/Domain/Statistics/Symbols.hs view
@@ -0,0 +1,133 @@+-----------------------------------------------------------------------------+-- Copyright 2020, Ideas project team. This file is distributed under the+-- terms of the Apache License 2.0. For more information, see the files+-- "LICENSE.txt" and "NOTICE.txt", which are included in the distribution.+-----------------------------------------------------------------------------+{- Implements symbols for the statistics domain. The only symbols that are+ - re-used are the true and the false symbols from the Logic1 module.+ -}+module Domain.Statistics.Symbols+ ( -- cChoiceSymbol+ cConstSymbol+ , cDataSymbol+ , componentSymbol+ -- , conclusionSymbol+ , derivedSymbol+ , falseSymbol+ , initialSymbol+ , leftSidedSymbol+ , rightSidedSymbol+ -- , sidedSymbol+ , singleSampleSymbol+ , tTestOneSymbol+ , tTestPairedSymbol+ , tTestTwoSymbol+ , rPearsonSymbol+ , anovaSymbol+ , chiSquaredSymbol+ , dependentSymbol, independentSymbol+ -- , testTypeSymbol+ , toOMSymbol+ , trueSymbol+ , twoSampleSymbol+ , twoSidedSymbol+ , zTestSymbol+ , rejectH0Symbol, acceptH0Symbol, dontRejectH0Symbol, rejectH1Symbol+ , acceptH1Symbol, dontRejectH1Symbol+ , toTermBool, fromTermBool+ ) where++import Ideas.Common.Rewriting.Term hiding (trueSymbol, falseSymbol)+import qualified Ideas.Text.OpenMath.Dictionary.Logic1 as L+import Ideas.Common.Id+import qualified Ideas.Text.OpenMath.Symbol as OM++componentSymbol :: Symbol+componentSymbol = newSymbol "stats.component"++-- cChoiceSymbol -- component data symbol+cConstSymbol -- component expression symbol, filtered in OM output+ , cDataSymbol -- component relation symbol, filtered in OM output+-- , conclusionSymbol -- All of these symbols up to intialSymbol are tied to a+ , derivedSymbol -- constructor from Domain.Statistics.ComponentSet+ , initialSymbol --+ , leftSidedSymbol --+ , rightSidedSymbol --+-- , sidedSymbol --+ , singleSampleSymbol --+ , tTestOneSymbol --+ , tTestPairedSymbol --+ , tTestTwoSymbol -- when the step is in the initial set+-- , testTypeSymbol -- when the step is derived+ , twoSampleSymbol+ , twoSidedSymbol+ , zTestSymbol+ , rPearsonSymbol+ , anovaSymbol+ , chiSquaredSymbol+ , dependentSymbol+ , independentSymbol :: Symbol++cDataSymbol = newSymbol "stats.cdata"+-- cChoiceSymbol = newSymbol "stats.cchoice"+cConstSymbol = newSymbol "stats.cconst"++--testTypeSymbol = newSymbol "stats.testtype"+--sidedSymbol = newSymbol "stats.sided"+--conclusionSymbol = newSymbol "stats.conclusion"++tTestOneSymbol = newSymbol "stats.ttestone"+tTestTwoSymbol = newSymbol "stats.ttesttwo"+tTestPairedSymbol = newSymbol "stats.ttestpaired"+zTestSymbol = newSymbol "stats.ztest"+rPearsonSymbol = newSymbol "stats.rpearson"+anovaSymbol = newSymbol "stats.anova"+chiSquaredSymbol = newSymbol "stats.chisquared"++twoSidedSymbol = newSymbol "stats.twosided"+leftSidedSymbol = newSymbol "stats.leftsided"+rightSidedSymbol = newSymbol "stats.rightsided"++singleSampleSymbol = newSymbol "stats.singlesample"+twoSampleSymbol = newSymbol "stats.twosample"++dependentSymbol = newSymbol "stats.dependent"+independentSymbol = newSymbol "stats.independent"++initialSymbol = newSymbol "stats.initial"+derivedSymbol = newSymbol "stats.derived"++trueSymbol, falseSymbol :: Symbol+trueSymbol = newSymbol L.trueSymbol+falseSymbol = newSymbol L.falseSymbol++toOMSymbol :: Symbol -> OM.Symbol+toOMSymbol s = idToSymbol (getId s)++rejectH0Symbol, acceptH0Symbol, dontRejectH0Symbol, rejectH1Symbol,+ acceptH1Symbol, dontRejectH1Symbol :: Symbol+rejectH0Symbol = newSymbol "stats.rejecth0" +acceptH0Symbol = newSymbol "stats.accepth0" +dontRejectH0Symbol = newSymbol "stats.dontrejecth0" +rejectH1Symbol = newSymbol "stats.rejecth1" +acceptH1Symbol = newSymbol "stats.accepth1" +dontRejectH1Symbol = newSymbol "stats.dontrejecth1"++-- Copied from Ideas.Encoding.OpenMathSuppoer+idToSymbol :: Id -> OM.Symbol+idToSymbol a+ | null (qualifiers a) =+ OM.extraSymbol (unqualified a)+ | otherwise =+ OM.makeSymbol (qualification a) (unqualified a)++toTermBool :: Bool -> Term+toTermBool b = symbol s+ where+ s = if b then trueSymbol else falseSymbol++fromTermBool :: Monad m => Term -> m Bool+fromTermBool (TCon s [])+ | s == trueSymbol = return True+ | s == falseSymbol = return False+fromTermBool _ = fail "Term is not a boolean"
+ src/Domain/Statistics/Views.hs view
@@ -0,0 +1,183 @@+-----------------------------------------------------------------------------+-- Copyright 2020, Ideas project team. This file is distributed under the+-- terms of the Apache License 2.0. For more information, see the files+-- "LICENSE.txt" and "NOTICE.txt", which are included in the distribution.+-----------------------------------------------------------------------------+module Domain.Statistics.Views+ ( Substitution, Substitutable(..), single+ , getSubstitution, substitute+ , checkSubstitution+ , substitutedView, evaluatedView+ ) where++import qualified Data.Map as M+import qualified Data.Set as S+import Data.List (intercalate)+import Data.Maybe+import Domain.Math.Data.Relation+import Domain.Math.Expr+import Domain.Math.Numeric.Views+import Domain.Statistics.ComponentSet+import Ideas.Common.View+import Ideas.Utils.Uniplate++-- idempotent, non-recursive substitutions+newtype Substitution = S (M.Map String Expr)++instance Show Substitution where+ show (S m) = "{" ++ intercalate ", " (map f (M.toList m)) ++ "}"+ where+ f (s, expr) = s ++ "=" ++ show expr++instance Semigroup Substitution where+ (<>) (S m) = flip (M.foldrWithKey extend) m++instance Monoid Substitution where+ mempty = S M.empty+ mappend = (<>)++member :: String -> Substitution -> Bool+member s (S m) = s `M.member` m++single :: String -> Expr -> Substitution+single s expr + | s `elem` vars expr = mempty+ | otherwise = S (M.singleton s expr)++extend :: String -> Expr -> Substitution -> Substitution+extend s expr sub@(S m)+ | s `elem` vars expr = mempty -- recursive expression+ | otherwise =+ case M.lookup s m of+ Nothing -> S (M.insert s (sub |-> expr) (M.map (single s expr |->) m))+ Just e + | e == (sub |-> expr) -> sub + -- inconsistent extension; should not happen.+ -- deleting s in m breaks monoid property for Substitution+ | otherwise -> S (M.delete s m) ++-- throws an error for non-idempotent substitutions+checkSubstitution :: Substitution -> Substitution +checkSubstitution sub@(S m)+ | disjoint (M.keysSet m) (vars (M.elems m)) = sub + | otherwise = error $ "Non-idempotent substitution: " ++ show sub+ where+ -- in Data.Set, since 0.5.11+ disjoint x y = S.null (S.intersection x y)++class Substitutable a where + (|->) :: Substitution -> a -> a+ vars :: a -> S.Set String++instance Substitutable a => Substitutable [a] where+ (|->) = map . (|->)+ vars = S.unions . map vars++instance Substitutable Expr where+ S xs |-> Var s = fromMaybe (Var s) (M.lookup s xs)+ sub |-> expr = descend (sub |->) expr+ + vars (Var s) = S.singleton s+ vars expr = S.unions (map vars (children expr))++instance Substitutable a => Substitutable (Relation a) where+ sub |-> r = fmap (sub |->) r+ vars r = vars (leftHandSide r) `S.union` vars (rightHandSide r)++instance Substitutable Component where+ sub |-> CExpr e = CExpr $ sub |-> e+ sub |-> CRelation r = CRelation $ sub |-> r+ _ |-> c = c+ + vars (CExpr e) = vars e+ vars (CRelation r) = vars r+ vars _ = S.empty++instance Substitutable ComponentSet where+ sub |-> cs = mapComponent (sub |->) cs+ vars = S.unions . map (vars . snd) . toList++getSubstitution :: ComponentSet -> Substitution+getSubstitution = addSigmaM . addSEM . make . resolveTestFormulaValue . addComponentVars+ where + -- if test value and test formula are both present (but differ because of + -- rounding), prefer the test value+ resolveTestFormulaValue cs + | cs `contains` TestValue && cs `contains` TestFormula = + delete TestFormula cs+ | otherwise = cs+ + make cs = mconcat (map f (xs ++ ys))+ where+ -- prefer initials in substitution+ xs = toList (initials cs)+ ys = filter notInital (toList (derived cs))+ notInital = (`notElem` map fst xs) . fst+ + f (cid, CExpr e) = single (getCId cid) e+ f (_, CRelation r) | relationType r == EqualTo =+ case leftHandSide r of+ Var v -> single v (rightHandSide r)+ _ -> mempty+ f _ = mempty+ + -- add formulas, but only if not already present in the substitution+ addSigmaM sub+ | "sigmaM" `member` sub = sub+ | otherwise = sub <> single "sigmaM" (Var "sigma" / sqrt (Var "n"))+ addSEM sub + | "SEM" `member` sub = sub+ | otherwise = sub <> single "SEM" (Var "s" / sqrt (Var "n"))++substitute :: ComponentSet -> ComponentSet+substitute cs = getSubstitution cs |-> cs++addComponentVars :: ComponentSet -> ComponentSet+addComponentVars cs = foldr (uncurry addComponentVar) (addPopulationMean cs) table+ where+ table = + [ (PopulationSdev, "sigma")+ , (SampleSdev, "s")+ , (SampleMean, "M")+ , (One SampleMean, "M1")+ , (Two SampleMean, "M2")+ , (One SampleSdev, "s1")+ , (Two SampleSdev, "s2")+ , (PopulationMean, "mu")+ ]++-- extract population mean from alternative hypothesis+addPopulationMean :: ComponentSet -> ComponentSet+addPopulationMean cs =+ case (get PopulationMean cs, getRelation AlternativeHypothesis cs) of+ (Nothing, Just rel) -> append PopulationMean (CExpr (rightHandSide rel)) cs+ _ -> cs++addComponentVar :: ComponentId -> String -> ComponentSet -> ComponentSet+addComponentVar cid var cs =+ case getRhsExpr cid cs of+ Just e | cs `doesNotContain` Other var ->+ -- to do: is it necessary to add (Other var) to the component set??+ append (Other var) (CExpr e) (single var e |-> cs)+ _ -> cs++----------------------------------------------------------+-- Substituted view++substitutedView :: View ComponentSet ComponentSet+substitutedView = makeView (Just . substitute) id++----------------------------------------------------------+-- Evaluated view++evaluatedView :: View ComponentSet ComponentSet+evaluatedView = makeView f id+ where+ f :: ComponentSet -> Maybe ComponentSet+ f cs = Just $ mapComponent evaluateComponent cs++evaluateComponent :: Component -> Component+evaluateComponent (CExpr x) = CExpr $ simplify doubleView x+evaluateComponent (CRelation r) = CRelation $ fmap (simplify doubleView) r+evaluateComponent c = c+
+ src/Main.hs view
@@ -0,0 +1,38 @@+-----------------------------------------------------------------------------+-- Copyright 2020, Ideas project team. This file is distributed under the+-- terms of the Apache License 2.0. For more information, see the files+-- "LICENSE.txt" and "NOTICE.txt", which are included in the distribution.+-----------------------------------------------------------------------------+module Main where++-- Include the Ideas framework+import Ideas.Common.Library+import Ideas.Main.Default++-- Include the hypothesis exercise and predicate service+import Domain.Hypothesis.Exercises+import Service.ConstraintServices++main :: IO ()+main = defaultMain ideasStatistics++ideasStatistics :: DomainReasoner+ideasStatistics =+ describe "Domain reasoner for statistics" $+ (newDomainReasoner "ideas.statistics")+ { exercises = myExercises+ , services = myServices+ }++myExercises :: [Some Exercise]+myExercises =+ [ -- hypothesis testing+ Some hypothesisExercise+ ]++myServices :: [Service]+myServices =+ concat [ metaServiceList ideasStatistics+ , filter ((/= newId "basic.diagnose") . getId) serviceList+ , constraintServiceList+ ]
+ src/Service/ConstraintServices.hs view
@@ -0,0 +1,93 @@+-----------------------------------------------------------------------------+-- Copyright 2020, Ideas project team. This file is distributed under the+-- terms of the Apache License 2.0. For more information, see the files+-- "LICENSE.txt" and "NOTICE.txt", which are included in the distribution.+-----------------------------------------------------------------------------+module Service.ConstraintServices+ ( constraintServiceList+ ) where++import Ideas.Common.Library+import Ideas.Service.Diagnose+import Ideas.Service.BasicServices+import Ideas.Service.State+import Ideas.Service.Types++-- Constraint services that are available+constraintServiceList :: [Service]+constraintServiceList = [diagnoseS, constraintsS, onehintS, allhintsS]++diagnoseS :: Service+diagnoseS = makeService "basic.diagnose"+ "Diagnose an expression submitted by a student. Possible diagnosis are \+ \Buggy (a common misconception was detected), NotEquivalent (something is \+ \wrong, but we don't know what), Similar (the expression is pretty similar \+ \to the last expression in the derivation), Expected (the submitted \+ \expression was anticipated by the strategy), Detour (the submitted \+ \expression was not expected by the strategy, but the applied rule was \+ \detected), and Correct (it is correct, but we don't know which rule was \+ \applied). Extended version for statistics domain: check predicates for NotEquiv." $+ diagnoseWithConstraints ::: tState .-> tContext .-> tMaybe tId .-> tDiagnosis++constraintsS :: Service+constraintsS = makeService "basic.constraints"+ "Check all constraints" $+ checkConstraints ::: tState .-> tList (tPair tConstraint tString)++onehintS :: Service+onehintS = makeService "basic.onehint"+ "Check constraints, and provide a hint" $+ onehint ::: tState .-> tError (Tag "elem" tHint)++allhintsS :: Service+allhintsS = makeService "basic.allhints"+ "Check constraints, and provide (multiple) hints" $+ allhints ::: tState .-> tError (tList tHint)++--------------------------------------------------------------------------------++diagnoseWithConstraints :: State a -> Context a -> Maybe Id -> Diagnosis a+diagnoseWithConstraints st ctx mid = f (diagnose st ctx mid)+ where+ f (NotEquivalent "") =+ case violations (exercise st) ctx of+ (n, msg):_ -> NotEquivalent $ violationMsg n msg+ [] | null (violations (exercise st) (stateContext st)) ->+ NotEquivalent []+ | otherwise -> -- special case: previous state is invalid (constraint violated)+ Correct (finished st) st+ f d = d++data Hint a = Violation (Constraint (Context a)) String+ | HintStep (StepInfo a) (State a)++tHint :: Type a (Hint a)+tHint = Iso (f <-> g) tp+ where+ tp = tPair tConstraint tString :|: tPair tStepInfo tState++ f (Left (c, msg)) = Violation c msg+ f (Right (stp, st)) = HintStep stp st+ g (Violation c msg) = Left (c, msg)+ g (HintStep stp st) = Right (stp, st)++allhints :: State a -> Either String [Hint a]+allhints st = + case violations (exercise st) (stateContext st) of + [] -> fmap (map (uncurry HintStep)) (allfirsts st)+ xs -> Right (map (uncurry Violation) xs)++onehint :: State a -> Either String (Hint a)+onehint st = + case allhints st of+ Left msg -> Left msg+ Right [] -> Left "no hint available"+ Right (h:_) -> Right h++checkConstraints :: State a -> [(Constraint (Context a), String)]+checkConstraints st = map f (constraints (exercise st))+ where+ f c = (c, show $ getResult c (stateContext st))++violationMsg :: Constraint a -> String -> String+violationMsg n msg = show n ++ if null msg then "" else ": " ++ msg