CARP-1.0/0000755000076400007640000000000011357714057012033 5ustar volmelnvolmelnCARP-1.0/GenPars.c0000664000076400007640000001114111356774573013546 0ustar volmelnvolmeln /* CARP: Clustering Algorithms' Referee Package Copyright (C) 2010 Volodymyr Melnykov and Ranjan Maitra This program is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version. This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details. You should have received a copy of the GNU General Public License along with this program. If not, see . Authors' contact information: Volodymyr Melnykov Ranjan Maitra Volodymyr.Melnykov@ndsu.edu maitra@iastate.edu Department of Statistics Department of Statistics North Dakota State University Iowa State University Fargo, ND 58102 Ames, IA 50011 */ #include #include #include #include "array.h" #include "overlap.h" /* genSigma : generates covariance matrix based on (p + 1) observations (unstable covariance matrix) * Parameters: * p - number of dimensions * VC - variance-covariance matrix */ void genSigma(int p, double **VC){ int i,j,k,n; double **x, *mu; n = p + 1; MAKE_MATRIX(x, n, p); MAKE_VECTOR(mu, p); anull(mu, p); for (i=0; i emax){ Anull(L, p, p); for (i=0; i= 1) | (PiLow <= 0)){ /* printf("Warning: PiLow is out of range... generated equal mixing proportions...\n"); */ for (k=0; k. Authors' contact information: Volodymyr Melnykov Ranjan Maitra Volodymyr.Melnykov@ndsu.edu maitra@iastate.edu Department of Statistics Department of Statistics North Dakota State University Iowa State University Fargo, ND 58102 Ames, IA 50011 */ #include #include #include #include "array.h" #include "overlap.h" /* genData: generates a dataset given mixture parameters * Parameters: * p - number of dimensions * K - number of components * Pi - vector of mixing proportions * Mu - set of mean vectors * S - set of covariance matrices * n - number of observations * Y - simulated dataset * Nk - cluster sizes */ void genData(int p, int K, double *Pi, double **Mu, double ***S, int n, double **Y, int *Nk){ int i, j, k; int curr; double **X; MAKE_MATRIX(X, n, p); rmulti(n, K, Pi, Nk); curr = 0; for (k=0; k. Authors' contact information: Volodymyr Melnykov Ranjan Maitra Volodymyr.Melnykov@ndsu.edu maitra@iastate.edu Department of Statistics Department of Statistics North Dakota State University Iowa State University Fargo, ND 58102 Ames, IA 50011 */ /* cephes_eigens.c * * This file was orignially called eigens.c and was found in the cephese * computational library. I edited it to bring it up to the ANSI C * standard and renamed it to cephes_eigens.c. * * --RR * July 2002 * * * Eigenvalues and eigenvectors of a real symmetric matrix * * * * SYNOPSIS: * * int n; * double A[n*(n+1)/2], EV[n*n], E[n]; * void cephes_eigens( A, EV, E, n ); * * * * DESCRIPTION: * * The algorithm is due to J. vonNeumann. * - - * A[] is a symmetric matrix stored in lower triangular form. * That is, A[ row, column ] = A[ (row*row+row)/2 + column ] * or equivalently with row and column interchanged. The * indices row and column run from 0 through n-1. * * EV[] is the output matrix of eigenvectors stored columnwise. * That is, the elements of each eigenvector appear in sequential * memory order. The jth element of the ith eigenvector is * EV[ n*i+j ] = EV[i][j]. * * E[] is the output matrix of eigenvalues. The ith element * of E corresponds to the ith eigenvector (the ith row of EV). * * On output, the matrix A will have been diagonalized and its * orginal contents are destroyed. * * ACCURACY: * * The error is controlled by an internal parameter called RANGE * which is set to 1e-10. After diagonalization, the * off-diagonal elements of A will have been reduced by * this factor. * * EEVOR MESSAGES: * * None. * */ /* Copyright 1973, 1991 by Stephen L. Moshier Copyleft version. */ #include "cephes_eigens.h" void cephes_eigens(double *A, double *EV, double *E, int N) { int IND, L, LL, LM, M, MM, MQ, I, J, IA, LQ; int IQ, IM, IL, NLI, NMI; double ANORM, ANORMX, AIA, THR, ALM, ALL, AMM, X, Y; double SINX, SINX2, COSX, COSX2, SINCS, AIL, AIM; double RLI, RMI; double sqrt(), fabs(); static double RANGE = 1.0e-10; /* 3.0517578e-5; */ /* Initialize identity matrix in EV[] */ for (J = 0; J < N * N; J++) EV[J] = 0.0; MM = 0; for (J = 0; J < N; J++) { EV[MM + J] = 1.0; MM += N; } ANORM = 0.0; for (I = 0; I < N; I++) { for (J = 0; J < N; J++) { if (I != J) { IA = I + (J * J + J) / 2; AIA = A[IA]; ANORM += AIA * AIA; } } } if (ANORM <= 0.0) goto done; ANORM = sqrt(ANORM + ANORM); ANORMX = ANORM * RANGE / N; THR = ANORM; while (THR > ANORMX) { THR = THR / N; do { /* while IND != 0 */ IND = 0; for (L = 0; L < N - 1; L++) { for (M = L + 1; M < N; M++) { MQ = (M * M + M) / 2; LM = L + MQ; ALM = A[LM]; if (fabs(ALM) < THR) continue; IND = 1; LQ = (L * L + L) / 2; LL = L + LQ; MM = M + MQ; ALL = A[LL]; AMM = A[MM]; X = (ALL - AMM) / 2.0; Y = -ALM / sqrt(ALM * ALM + X * X); if (X < 0.0) Y = -Y; SINX = Y/sqrt(2.0*(1.0+sqrt(1.0-Y*Y))); SINX2 = SINX * SINX; COSX = sqrt(1.0 - SINX2); COSX2 = COSX * COSX; SINCS = SINX * COSX; /* ROTATE L AND M COLUMNS */ for (I = 0; I < N; I++) { IQ = (I * I + I) / 2; if ((I != M) && (I != L)) { if (I > M) IM = M + IQ; else IM = I + MQ; if (I >= L) IL = L + IQ; else IL = I + LQ; AIL = A[IL]; AIM = A[IM]; X = AIL * COSX - AIM * SINX; A[IM] = AIL * SINX + AIM * COSX; A[IL] = X; } NLI = N * L + I; NMI = N * M + I; RLI = EV[NLI]; RMI = EV[NMI]; EV[NLI] = RLI * COSX - RMI * SINX; EV[NMI] = RLI * SINX + RMI * COSX; } X = 2.0 * ALM * SINCS; A[LL] = ALL * COSX2 + AMM * SINX2 - X; A[MM] = ALL * SINX2 + AMM * COSX2 + X; A[LM] = (ALL - AMM) * SINCS + ALM * (COSX2 - SINX2); } /* for M=L+1 to N-1 */ } /* for L=0 to N-2 */ } while (IND != 0); } /* while THR > ANORMX */ done: /* Extract eigenvalues from the reduced matrix */ L = 0; for (J = 1; J <= N; J++) { L = L + J; E[J - 1] = A[L - 1]; } } CARP-1.0/COPYING0000664000076400007640000010451311356252646013074 0ustar volmelnvolmeln GNU GENERAL PUBLIC LICENSE Version 3, 29 June 2007 Copyright (C) 2007 Free Software Foundation, Inc. Everyone is permitted to copy and distribute verbatim copies of this license document, but changing it is not allowed. Preamble The GNU General Public License is a free, copyleft license for software and other kinds of works. The licenses for most software and other practical works are designed to take away your freedom to share and change the works. By contrast, the GNU General Public License is intended to guarantee your freedom to share and change all versions of a program--to make sure it remains free software for all its users. We, the Free Software Foundation, use the GNU General Public License for most of our software; it applies also to any other work released this way by its authors. You can apply it to your programs, too. When we speak of free software, we are referring to freedom, not price. Our General Public Licenses are designed to make sure that you have the freedom to distribute copies of free software (and charge for them if you wish), that you receive source code or can get it if you want it, that you can change the software or use pieces of it in new free programs, and that you know you can do these things. To protect your rights, we need to prevent others from denying you these rights or asking you to surrender the rights. Therefore, you have certain responsibilities if you distribute copies of the software, or if you modify it: responsibilities to respect the freedom of others. For example, if you distribute copies of such a program, whether gratis or for a fee, you must pass on to the recipients the same freedoms that you received. You must make sure that they, too, receive or can get the source code. And you must show them these terms so they know their rights. Developers that use the GNU GPL protect your rights with two steps: (1) assert copyright on the software, and (2) offer you this License giving you legal permission to copy, distribute and/or modify it. For the developers' and authors' protection, the GPL clearly explains that there is no warranty for this free software. For both users' and authors' sake, the GPL requires that modified versions be marked as changed, so that their problems will not be attributed erroneously to authors of previous versions. Some devices are designed to deny users access to install or run modified versions of the software inside them, although the manufacturer can do so. This is fundamentally incompatible with the aim of protecting users' freedom to change the software. The systematic pattern of such abuse occurs in the area of products for individuals to use, which is precisely where it is most unacceptable. Therefore, we have designed this version of the GPL to prohibit the practice for those products. If such problems arise substantially in other domains, we stand ready to extend this provision to those domains in future versions of the GPL, as needed to protect the freedom of users. Finally, every program is threatened constantly by software patents. States should not allow patents to restrict development and use of software on general-purpose computers, but in those that do, we wish to avoid the special danger that patents applied to a free program could make it effectively proprietary. To prevent this, the GPL assures that patents cannot be used to render the program non-free. The precise terms and conditions for copying, distribution and modification follow. TERMS AND CONDITIONS 0. Definitions. "This License" refers to version 3 of the GNU General Public License. "Copyright" also means copyright-like laws that apply to other kinds of works, such as semiconductor masks. "The Program" refers to any copyrightable work licensed under this License. Each licensee is addressed as "you". "Licensees" and "recipients" may be individuals or organizations. To "modify" a work means to copy from or adapt all or part of the work in a fashion requiring copyright permission, other than the making of an exact copy. The resulting work is called a "modified version" of the earlier work or a work "based on" the earlier work. A "covered work" means either the unmodified Program or a work based on the Program. To "propagate" a work means to do anything with it that, without permission, would make you directly or secondarily liable for infringement under applicable copyright law, except executing it on a computer or modifying a private copy. Propagation includes copying, distribution (with or without modification), making available to the public, and in some countries other activities as well. To "convey" a work means any kind of propagation that enables other parties to make or receive copies. Mere interaction with a user through a computer network, with no transfer of a copy, is not conveying. An interactive user interface displays "Appropriate Legal Notices" to the extent that it includes a convenient and prominently visible feature that (1) displays an appropriate copyright notice, and (2) tells the user that there is no warranty for the work (except to the extent that warranties are provided), that licensees may convey the work under this License, and how to view a copy of this License. If the interface presents a list of user commands or options, such as a menu, a prominent item in the list meets this criterion. 1. Source Code. The "source code" for a work means the preferred form of the work for making modifications to it. "Object code" means any non-source form of a work. A "Standard Interface" means an interface that either is an official standard defined by a recognized standards body, or, in the case of interfaces specified for a particular programming language, one that is widely used among developers working in that language. The "System Libraries" of an executable work include anything, other than the work as a whole, that (a) is included in the normal form of packaging a Major Component, but which is not part of that Major Component, and (b) serves only to enable use of the work with that Major Component, or to implement a Standard Interface for which an implementation is available to the public in source code form. A "Major Component", in this context, means a major essential component (kernel, window system, and so on) of the specific operating system (if any) on which the executable work runs, or a compiler used to produce the work, or an object code interpreter used to run it. The "Corresponding Source" for a work in object code form means all the source code needed to generate, install, and (for an executable work) run the object code and to modify the work, including scripts to control those activities. However, it does not include the work's System Libraries, or general-purpose tools or generally available free programs which are used unmodified in performing those activities but which are not part of the work. For example, Corresponding Source includes interface definition files associated with source files for the work, and the source code for shared libraries and dynamically linked subprograms that the work is specifically designed to require, such as by intimate data communication or control flow between those subprograms and other parts of the work. The Corresponding Source need not include anything that users can regenerate automatically from other parts of the Corresponding Source. The Corresponding Source for a work in source code form is that same work. 2. Basic Permissions. All rights granted under this License are granted for the term of copyright on the Program, and are irrevocable provided the stated conditions are met. This License explicitly affirms your unlimited permission to run the unmodified Program. The output from running a covered work is covered by this License only if the output, given its content, constitutes a covered work. This License acknowledges your rights of fair use or other equivalent, as provided by copyright law. You may make, run and propagate covered works that you do not convey, without conditions so long as your license otherwise remains in force. You may convey covered works to others for the sole purpose of having them make modifications exclusively for you, or provide you with facilities for running those works, provided that you comply with the terms of this License in conveying all material for which you do not control copyright. Those thus making or running the covered works for you must do so exclusively on your behalf, under your direction and control, on terms that prohibit them from making any copies of your copyrighted material outside their relationship with you. Conveying under any other circumstances is permitted solely under the conditions stated below. Sublicensing is not allowed; section 10 makes it unnecessary. 3. Protecting Users' Legal Rights From Anti-Circumvention Law. No covered work shall be deemed part of an effective technological measure under any applicable law fulfilling obligations under article 11 of the WIPO copyright treaty adopted on 20 December 1996, or similar laws prohibiting or restricting circumvention of such measures. When you convey a covered work, you waive any legal power to forbid circumvention of technological measures to the extent such circumvention is effected by exercising rights under this License with respect to the covered work, and you disclaim any intention to limit operation or modification of the work as a means of enforcing, against the work's users, your or third parties' legal rights to forbid circumvention of technological measures. 4. Conveying Verbatim Copies. You may convey verbatim copies of the Program's source code as you receive it, in any medium, provided that you conspicuously and appropriately publish on each copy an appropriate copyright notice; keep intact all notices stating that this License and any non-permissive terms added in accord with section 7 apply to the code; keep intact all notices of the absence of any warranty; and give all recipients a copy of this License along with the Program. You may charge any price or no price for each copy that you convey, and you may offer support or warranty protection for a fee. 5. Conveying Modified Source Versions. You may convey a work based on the Program, or the modifications to produce it from the Program, in the form of source code under the terms of section 4, provided that you also meet all of these conditions: a) The work must carry prominent notices stating that you modified it, and giving a relevant date. b) The work must carry prominent notices stating that it is released under this License and any conditions added under section 7. This requirement modifies the requirement in section 4 to "keep intact all notices". c) You must license the entire work, as a whole, under this License to anyone who comes into possession of a copy. This License will therefore apply, along with any applicable section 7 additional terms, to the whole of the work, and all its parts, regardless of how they are packaged. This License gives no permission to license the work in any other way, but it does not invalidate such permission if you have separately received it. d) If the work has interactive user interfaces, each must display Appropriate Legal Notices; however, if the Program has interactive interfaces that do not display Appropriate Legal Notices, your work need not make them do so. A compilation of a covered work with other separate and independent works, which are not by their nature extensions of the covered work, and which are not combined with it such as to form a larger program, in or on a volume of a storage or distribution medium, is called an "aggregate" if the compilation and its resulting copyright are not used to limit the access or legal rights of the compilation's users beyond what the individual works permit. Inclusion of a covered work in an aggregate does not cause this License to apply to the other parts of the aggregate. 6. Conveying Non-Source Forms. You may convey a covered work in object code form under the terms of sections 4 and 5, provided that you also convey the machine-readable Corresponding Source under the terms of this License, in one of these ways: a) Convey the object code in, or embodied in, a physical product (including a physical distribution medium), accompanied by the Corresponding Source fixed on a durable physical medium customarily used for software interchange. b) Convey the object code in, or embodied in, a physical product (including a physical distribution medium), accompanied by a written offer, valid for at least three years and valid for as long as you offer spare parts or customer support for that product model, to give anyone who possesses the object code either (1) a copy of the Corresponding Source for all the software in the product that is covered by this License, on a durable physical medium customarily used for software interchange, for a price no more than your reasonable cost of physically performing this conveying of source, or (2) access to copy the Corresponding Source from a network server at no charge. c) Convey individual copies of the object code with a copy of the written offer to provide the Corresponding Source. This alternative is allowed only occasionally and noncommercially, and only if you received the object code with such an offer, in accord with subsection 6b. d) Convey the object code by offering access from a designated place (gratis or for a charge), and offer equivalent access to the Corresponding Source in the same way through the same place at no further charge. You need not require recipients to copy the Corresponding Source along with the object code. If the place to copy the object code is a network server, the Corresponding Source may be on a different server (operated by you or a third party) that supports equivalent copying facilities, provided you maintain clear directions next to the object code saying where to find the Corresponding Source. Regardless of what server hosts the Corresponding Source, you remain obligated to ensure that it is available for as long as needed to satisfy these requirements. e) Convey the object code using peer-to-peer transmission, provided you inform other peers where the object code and Corresponding Source of the work are being offered to the general public at no charge under subsection 6d. A separable portion of the object code, whose source code is excluded from the Corresponding Source as a System Library, need not be included in conveying the object code work. A "User Product" is either (1) a "consumer product", which means any tangible personal property which is normally used for personal, family, or household purposes, or (2) anything designed or sold for incorporation into a dwelling. In determining whether a product is a consumer product, doubtful cases shall be resolved in favor of coverage. For a particular product received by a particular user, "normally used" refers to a typical or common use of that class of product, regardless of the status of the particular user or of the way in which the particular user actually uses, or expects or is expected to use, the product. A product is a consumer product regardless of whether the product has substantial commercial, industrial or non-consumer uses, unless such uses represent the only significant mode of use of the product. "Installation Information" for a User Product means any methods, procedures, authorization keys, or other information required to install and execute modified versions of a covered work in that User Product from a modified version of its Corresponding Source. The information must suffice to ensure that the continued functioning of the modified object code is in no case prevented or interfered with solely because modification has been made. If you convey an object code work under this section in, or with, or specifically for use in, a User Product, and the conveying occurs as part of a transaction in which the right of possession and use of the User Product is transferred to the recipient in perpetuity or for a fixed term (regardless of how the transaction is characterized), the Corresponding Source conveyed under this section must be accompanied by the Installation Information. But this requirement does not apply if neither you nor any third party retains the ability to install modified object code on the User Product (for example, the work has been installed in ROM). The requirement to provide Installation Information does not include a requirement to continue to provide support service, warranty, or updates for a work that has been modified or installed by the recipient, or for the User Product in which it has been modified or installed. Access to a network may be denied when the modification itself materially and adversely affects the operation of the network or violates the rules and protocols for communication across the network. Corresponding Source conveyed, and Installation Information provided, in accord with this section must be in a format that is publicly documented (and with an implementation available to the public in source code form), and must require no special password or key for unpacking, reading or copying. 7. Additional Terms. "Additional permissions" are terms that supplement the terms of this License by making exceptions from one or more of its conditions. Additional permissions that are applicable to the entire Program shall be treated as though they were included in this License, to the extent that they are valid under applicable law. If additional permissions apply only to part of the Program, that part may be used separately under those permissions, but the entire Program remains governed by this License without regard to the additional permissions. When you convey a copy of a covered work, you may at your option remove any additional permissions from that copy, or from any part of it. (Additional permissions may be written to require their own removal in certain cases when you modify the work.) You may place additional permissions on material, added by you to a covered work, for which you have or can give appropriate copyright permission. Notwithstanding any other provision of this License, for material you add to a covered work, you may (if authorized by the copyright holders of that material) supplement the terms of this License with terms: a) Disclaiming warranty or limiting liability differently from the terms of sections 15 and 16 of this License; or b) Requiring preservation of specified reasonable legal notices or author attributions in that material or in the Appropriate Legal Notices displayed by works containing it; or c) Prohibiting misrepresentation of the origin of that material, or requiring that modified versions of such material be marked in reasonable ways as different from the original version; or d) Limiting the use for publicity purposes of names of licensors or authors of the material; or e) Declining to grant rights under trademark law for use of some trade names, trademarks, or service marks; or f) Requiring indemnification of licensors and authors of that material by anyone who conveys the material (or modified versions of it) with contractual assumptions of liability to the recipient, for any liability that these contractual assumptions directly impose on those licensors and authors. All other non-permissive additional terms are considered "further restrictions" within the meaning of section 10. If the Program as you received it, or any part of it, contains a notice stating that it is governed by this License along with a term that is a further restriction, you may remove that term. If a license document contains a further restriction but permits relicensing or conveying under this License, you may add to a covered work material governed by the terms of that license document, provided that the further restriction does not survive such relicensing or conveying. If you add terms to a covered work in accord with this section, you must place, in the relevant source files, a statement of the additional terms that apply to those files, or a notice indicating where to find the applicable terms. Additional terms, permissive or non-permissive, may be stated in the form of a separately written license, or stated as exceptions; the above requirements apply either way. 8. Termination. You may not propagate or modify a covered work except as expressly provided under this License. Any attempt otherwise to propagate or modify it is void, and will automatically terminate your rights under this License (including any patent licenses granted under the third paragraph of section 11). However, if you cease all violation of this License, then your license from a particular copyright holder is reinstated (a) provisionally, unless and until the copyright holder explicitly and finally terminates your license, and (b) permanently, if the copyright holder fails to notify you of the violation by some reasonable means prior to 60 days after the cessation. Moreover, your license from a particular copyright holder is reinstated permanently if the copyright holder notifies you of the violation by some reasonable means, this is the first time you have received notice of violation of this License (for any work) from that copyright holder, and you cure the violation prior to 30 days after your receipt of the notice. Termination of your rights under this section does not terminate the licenses of parties who have received copies or rights from you under this License. If your rights have been terminated and not permanently reinstated, you do not qualify to receive new licenses for the same material under section 10. 9. Acceptance Not Required for Having Copies. You are not required to accept this License in order to receive or run a copy of the Program. Ancillary propagation of a covered work occurring solely as a consequence of using peer-to-peer transmission to receive a copy likewise does not require acceptance. However, nothing other than this License grants you permission to propagate or modify any covered work. These actions infringe copyright if you do not accept this License. Therefore, by modifying or propagating a covered work, you indicate your acceptance of this License to do so. 10. Automatic Licensing of Downstream Recipients. Each time you convey a covered work, the recipient automatically receives a license from the original licensors, to run, modify and propagate that work, subject to this License. You are not responsible for enforcing compliance by third parties with this License. An "entity transaction" is a transaction transferring control of an organization, or substantially all assets of one, or subdividing an organization, or merging organizations. If propagation of a covered work results from an entity transaction, each party to that transaction who receives a copy of the work also receives whatever licenses to the work the party's predecessor in interest had or could give under the previous paragraph, plus a right to possession of the Corresponding Source of the work from the predecessor in interest, if the predecessor has it or can get it with reasonable efforts. You may not impose any further restrictions on the exercise of the rights granted or affirmed under this License. For example, you may not impose a license fee, royalty, or other charge for exercise of rights granted under this License, and you may not initiate litigation (including a cross-claim or counterclaim in a lawsuit) alleging that any patent claim is infringed by making, using, selling, offering for sale, or importing the Program or any portion of it. 11. Patents. A "contributor" is a copyright holder who authorizes use under this License of the Program or a work on which the Program is based. The work thus licensed is called the contributor's "contributor version". A contributor's "essential patent claims" are all patent claims owned or controlled by the contributor, whether already acquired or hereafter acquired, that would be infringed by some manner, permitted by this License, of making, using, or selling its contributor version, but do not include claims that would be infringed only as a consequence of further modification of the contributor version. For purposes of this definition, "control" includes the right to grant patent sublicenses in a manner consistent with the requirements of this License. Each contributor grants you a non-exclusive, worldwide, royalty-free patent license under the contributor's essential patent claims, to make, use, sell, offer for sale, import and otherwise run, modify and propagate the contents of its contributor version. In the following three paragraphs, a "patent license" is any express agreement or commitment, however denominated, not to enforce a patent (such as an express permission to practice a patent or covenant not to sue for patent infringement). To "grant" such a patent license to a party means to make such an agreement or commitment not to enforce a patent against the party. If you convey a covered work, knowingly relying on a patent license, and the Corresponding Source of the work is not available for anyone to copy, free of charge and under the terms of this License, through a publicly available network server or other readily accessible means, then you must either (1) cause the Corresponding Source to be so available, or (2) arrange to deprive yourself of the benefit of the patent license for this particular work, or (3) arrange, in a manner consistent with the requirements of this License, to extend the patent license to downstream recipients. "Knowingly relying" means you have actual knowledge that, but for the patent license, your conveying the covered work in a country, or your recipient's use of the covered work in a country, would infringe one or more identifiable patents in that country that you have reason to believe are valid. If, pursuant to or in connection with a single transaction or arrangement, you convey, or propagate by procuring conveyance of, a covered work, and grant a patent license to some of the parties receiving the covered work authorizing them to use, propagate, modify or convey a specific copy of the covered work, then the patent license you grant is automatically extended to all recipients of the covered work and works based on it. A patent license is "discriminatory" if it does not include within the scope of its coverage, prohibits the exercise of, or is conditioned on the non-exercise of one or more of the rights that are specifically granted under this License. You may not convey a covered work if you are a party to an arrangement with a third party that is in the business of distributing software, under which you make payment to the third party based on the extent of your activity of conveying the work, and under which the third party grants, to any of the parties who would receive the covered work from you, a discriminatory patent license (a) in connection with copies of the covered work conveyed by you (or copies made from those copies), or (b) primarily for and in connection with specific products or compilations that contain the covered work, unless you entered into that arrangement, or that patent license was granted, prior to 28 March 2007. Nothing in this License shall be construed as excluding or limiting any implied license or other defenses to infringement that may otherwise be available to you under applicable patent law. 12. No Surrender of Others' Freedom. If conditions are imposed on you (whether by court order, agreement or otherwise) that contradict the conditions of this License, they do not excuse you from the conditions of this License. If you cannot convey a covered work so as to satisfy simultaneously your obligations under this License and any other pertinent obligations, then as a consequence you may not convey it at all. For example, if you agree to terms that obligate you to collect a royalty for further conveying from those to whom you convey the Program, the only way you could satisfy both those terms and this License would be to refrain entirely from conveying the Program. 13. Use with the GNU Affero General Public License. Notwithstanding any other provision of this License, you have permission to link or combine any covered work with a work licensed under version 3 of the GNU Affero General Public License into a single combined work, and to convey the resulting work. The terms of this License will continue to apply to the part which is the covered work, but the special requirements of the GNU Affero General Public License, section 13, concerning interaction through a network will apply to the combination as such. 14. Revised Versions of this License. The Free Software Foundation may publish revised and/or new versions of the GNU General Public License from time to time. Such new versions will be similar in spirit to the present version, but may differ in detail to address new problems or concerns. Each version is given a distinguishing version number. If the Program specifies that a certain numbered version of the GNU General Public License "or any later version" applies to it, you have the option of following the terms and conditions either of that numbered version or of any later version published by the Free Software Foundation. If the Program does not specify a version number of the GNU General Public License, you may choose any version ever published by the Free Software Foundation. If the Program specifies that a proxy can decide which future versions of the GNU General Public License can be used, that proxy's public statement of acceptance of a version permanently authorizes you to choose that version for the Program. Later license versions may give you additional or different permissions. However, no additional obligations are imposed on any author or copyright holder as a result of your choosing to follow a later version. 15. Disclaimer of Warranty. THERE IS NO WARRANTY FOR THE PROGRAM, TO THE EXTENT PERMITTED BY APPLICABLE LAW. EXCEPT WHEN OTHERWISE STATED IN WRITING THE COPYRIGHT HOLDERS AND/OR OTHER PARTIES PROVIDE THE PROGRAM "AS IS" WITHOUT WARRANTY OF ANY KIND, EITHER EXPRESSED OR IMPLIED, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE. THE ENTIRE RISK AS TO THE QUALITY AND PERFORMANCE OF THE PROGRAM IS WITH YOU. SHOULD THE PROGRAM PROVE DEFECTIVE, YOU ASSUME THE COST OF ALL NECESSARY SERVICING, REPAIR OR CORRECTION. 16. Limitation of Liability. IN NO EVENT UNLESS REQUIRED BY APPLICABLE LAW OR AGREED TO IN WRITING WILL ANY COPYRIGHT HOLDER, OR ANY OTHER PARTY WHO MODIFIES AND/OR CONVEYS THE PROGRAM AS PERMITTED ABOVE, BE LIABLE TO YOU FOR DAMAGES, INCLUDING ANY GENERAL, SPECIAL, INCIDENTAL OR CONSEQUENTIAL DAMAGES ARISING OUT OF THE USE OR INABILITY TO USE THE PROGRAM (INCLUDING BUT NOT LIMITED TO LOSS OF DATA OR DATA BEING RENDERED INACCURATE OR LOSSES SUSTAINED BY YOU OR THIRD PARTIES OR A FAILURE OF THE PROGRAM TO OPERATE WITH ANY OTHER PROGRAMS), EVEN IF SUCH HOLDER OR OTHER PARTY HAS BEEN ADVISED OF THE POSSIBILITY OF SUCH DAMAGES. 17. Interpretation of Sections 15 and 16. If the disclaimer of warranty and limitation of liability provided above cannot be given local legal effect according to their terms, reviewing courts shall apply local law that most closely approximates an absolute waiver of all civil liability in connection with the Program, unless a warranty or assumption of liability accompanies a copy of the Program in return for a fee. END OF TERMS AND CONDITIONS How to Apply These Terms to Your New Programs If you develop a new program, and you want it to be of the greatest possible use to the public, the best way to achieve this is to make it free software which everyone can redistribute and change under these terms. To do so, attach the following notices to the program. It is safest to attach them to the start of each source file to most effectively state the exclusion of warranty; and each file should have at least the "copyright" line and a pointer to where the full notice is found. Copyright (C) This program is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version. This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details. You should have received a copy of the GNU General Public License along with this program. If not, see . Also add information on how to contact you by electronic and paper mail. If the program does terminal interaction, make it output a short notice like this when it starts in an interactive mode: Copyright (C) This program comes with ABSOLUTELY NO WARRANTY; for details type `show w'. This is free software, and you are welcome to redistribute it under certain conditions; type `show c' for details. The hypothetical commands `show w' and `show c' should show the appropriate parts of the General Public License. Of course, your program's commands might be different; for a GUI interface, you would use an "about box". You should also get your employer (if you work as a programmer) or school, if any, to sign a "copyright disclaimer" for the program, if necessary. For more information on this, and how to apply and follow the GNU GPL, see . The GNU General Public License does not permit incorporating your program into proprietary programs. If your program is a subroutine library, you may consider it more useful to permit linking proprietary applications with the library. If this is what you want to do, use the GNU Lesser General Public License instead of this License. But first, please read . CARP-1.0/symmEigValDec.c0000664000076400007640000000456111356777403014703 0ustar volmelnvolmeln /* CARP: Clustering Algorithms' Referee Package Copyright (C) 2010 Volodymyr Melnykov and Ranjan Maitra This program is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version. This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details. You should have received a copy of the GNU General Public License along with this program. If not, see . Authors' contact information: Volodymyr Melnykov Ranjan Maitra Volodymyr.Melnykov@ndsu.edu maitra@iastate.edu Department of Statistics Department of Statistics North Dakota State University Iowa State University Fargo, ND 58102 Ames, IA 50011 */ #include #include #include #include "array.h" #include "cephes_eigens.h" void cephes_symmeigens_down(int p, double *eval, double **A, double (*determinant)) /* input: p - dimension of matrix A - pxp symmetric matrix (only lower triangle used) destroyed on return returns: eval - vector of eigenvalues (in ascending order) evec - matrix of eigenvectors determinant - determinant of the symmetric matrix (as calculated by the product of the eigenvalues */ { int i, j; double *As, *Evec, *Evalues; MAKE_VECTOR(As, p * (p + 1) / 2); for (i = 0; i < p; i++) { for (j = 0; j <= i; j++) As[(i * i + i)/2 + j] = A[i][j]; } MAKE_VECTOR(Evec, p * p); MAKE_VECTOR(Evalues, p); cephes_eigens(As, Evec, Evalues, p); for (i = 0; i < p; i++){ eval[i] = Evalues[p - i - 1]; } for (i = 0; i < p; i++) { for (j = 0; j < p; j++){ A[j][p-i-1] = Evec[p * i + j]; } } /* printf("EigenValues:\n"); for (i = 0; i < p; i++){ printf("EigVal: %lf \n", eval[i]); } printf("EigenVectors:\n"); for (i = 0; i < p; i++) { for (j = 0; j < p; j++){ printf("%lf ", A[i][j]); } printf("\n"); } */ (*determinant)=1.0; for (i = 0; i < p; i++) (*determinant)*=eval[i]; FREE_VECTOR(As); FREE_VECTOR(Evalues); FREE_VECTOR(Evec); return; } CARP-1.0/unittest.c0000664000076400007640000001176311356777421014073 0ustar volmelnvolmeln /* CARP: Clustering Algorithms' Referee Package Copyright (C) 2010 Volodymyr Melnykov and Ranjan Maitra This program is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version. This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details. You should have received a copy of the GNU General Public License along with this program. If not, see . Authors' contact information: Volodymyr Melnykov Ranjan Maitra Volodymyr.Melnykov@ndsu.edu maitra@iastate.edu Department of Statistics Department of Statistics North Dakota State University Iowa State University Fargo, ND 58102 Ames, IA 50011 */ #include #include #include #include #include #include #include "array.h" #include "overlap.h" #include "minunit.h" int getopt(int argc, char *const argv[], const char *optstring); extern char *optarg; extern int optind, optopt; int main(int argc, char *argv[]) { int k, sum; int p, K, resN, method, lim, fail, sph, g, nf; double emax, PiLow, Ubound; double BarOmega, MaxOmega; int *rcMax, *Nk; double *Pi, *pars; double **Mu, **OmegaMap, **Y; double ***S; double tol; int cond, cond1, cond2; MAKE_VECTOR(pars, 2); /* DEFAULT PARAMETERS */ p = 5; K = 7; method = 0; /* 0 - average; 1 - maximum */ BarOmega = 0.0; MaxOmega = 0.0; sph = 0; /* 0 - nonspherical, 1 - spherical */ emax = 0.90; /* Max eccentricity */ PiLow = 1.0; /* Smallest mixing proportion */ Ubound = 1.0; /* Means are distributed by U(0,Ubound) */ resN = 100; /* Max # of dataset resimulations */ pars[0] = 0.0000001; /* eps */ pars[1] = 0.0000001; /* acc (Davies) */ lim = 1000000; /* lim (Davies) */ g = 100; /* sample size for generated datasets */ nf = 1; /* # of simulated mixtures */ char PIfname[100] = "TEST/Pi.test"; char MUfname[100] = "TEST/Mu.test"; char Sfname[100] = "TEST/LTSigma.test"; srand(time(NULL)); MAKE_VECTOR(Pi, K); MAKE_MATRIX(Mu, K, p); MAKE_3ARRAY(S, K, p, p); MAKE_MATRIX(OmegaMap, K, K); MAKE_VECTOR(rcMax, 2); MAKE_MATRIX(Y, g, p); MAKE_VECTOR(Nk, K); /* RUNNING A CORRESPONDING FUNCTION */ tol = 0.000001; /* FUNCTION freadParameters */ freadParameters(p, K, Pi, Mu, S, PIfname, MUfname, Sfname, 1); cond = (p == 5) && (K == 7) && (fabs(Pi[0] - 0.1831131) < tol) && (fabs(Pi[6] - 0.1569579) < tol) && (fabs(Mu[0][0] - 0.0012522) < tol) && (fabs(Mu[1][1] - 0.6093170) < tol) && (fabs(Mu[2][2] - 0.3875023) < tol) && (fabs(Mu[3][3] - 0.9270482) < tol) && (fabs(Mu[4][4] - 0.7657742) < tol) && (fabs(Mu[6][4] - 0.7767524) < tol) && (fabs(S[0][0][0] - 0.0162117) < tol) && (fabs(S[1][0][1] + 0.0005859) < tol) && (fabs(S[2][0][2] + 0.0060184) < tol) && (fabs(S[4][3][0] - 0.0047450) < tol) && (fabs(S[5][0][4] - 0.0008476) < tol) && (fabs(S[6][4][4] - 0.0208782) < tol); printf("function \"freadParameters\" - "); func_test(cond, "OK\n", "Test failed\n"); /* FUNCTION ExactOverlap */ ExactOverlap(p, K, Pi, Mu, S, pars, lim, OmegaMap, &BarOmega, &MaxOmega, rcMax); cond = (fabs(BarOmega - 0.017414) < tol) && (fabs(MaxOmega - 0.198721) < tol); printf("function \"ExactOverlap\" - "); func_test(cond, "OK\n", "Test failed\n"); /* FUNCTION genData */ Anull(Y, g, p); genData(p, K, Pi, Mu, S, g, Y, Nk); sum = 0; for (k=0; k < K; k++){ sum = sum + Nk[k]; } cond = (g == sum) && (Y[0][0] != 0.0) && (Y[g-1][p-1] != 0.0); printf("function \"genData\" - "); func_test(cond, "OK\n", "Test failed\n"); /* FUNCTION OmegaClust */ // test 1 sph = 0; method = 1; OmegaClust(0.15, method, p, K, PiLow, Ubound, emax, pars, lim, resN, sph, Pi, Mu, S, OmegaMap, &BarOmega, &MaxOmega, rcMax, &fail); cond1 = (fabs(0.15 - MaxOmega) < tol); // test 2 PiLow = 0.1; sph = 1; method = 0; OmegaClust(0.03, method, p, K, PiLow, Ubound, emax, pars, lim, resN, sph, Pi, Mu, S, OmegaMap, &BarOmega, &MaxOmega, rcMax, &fail); cond2 = (fabs(0.03 - BarOmega) < tol) && (fabs(S[0][1][0]) < tol) && (fabs(S[1][2][1]) < tol); cond = (cond1 && cond2); printf("function \"OmegaClust\" - "); func_test(cond, "OK\n", "Test failed\n"); /* FUNCTION OmegaClust */ sph = 0; BarOmega = 0.03; MaxOmega = 0.15; OmegaBarOmegaMax(p, K, PiLow, Ubound, emax, pars, lim, resN, sph, Pi, Mu, S, OmegaMap, &BarOmega, &MaxOmega, rcMax, &fail); cond = (fabs(BarOmega - 0.03) < tol) && (fabs(MaxOmega - 0.15) < tol); printf("function \"OmegaBarOmegaMax\" - "); func_test(cond, "OK\n", "Test failed\n"); FREE_MATRIX(Y); FREE_VECTOR(Nk); FREE_VECTOR(Pi); FREE_MATRIX(Mu); FREE_3ARRAY(S); FREE_MATRIX(OmegaMap); FREE_VECTOR(rcMax); FREE_VECTOR(pars); return 0; } CARP-1.0/overlap.h0000644000076400007640000001045011356775076013663 0ustar volmelnvolmeln /* CARP: Clustering Algorithms' Referee Package Copyright (C) 2010 Volodymyr Melnykov and Ranjan Maitra This program is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version. This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details. You should have received a copy of the GNU General Public License along with this program. If not, see . Authors' contact information: Volodymyr Melnykov Ranjan Maitra Volodymyr.Melnykov@ndsu.edu maitra@iastate.edu Department of Statistics Department of Statistics North Dakota State University Iowa State University Fargo, ND 58102 Ames, IA 50011 */ #ifndef OVERLAP_H #define OVERLAP_H void ExactOverlap(int p, int K, double *Pi, double **Mu, double ***S, double *pars, int lim, double **OmegaMap, double (*BarOmega), double (*MaxOmega), int *rcMax); void OmegaClust(double Omega, int method, int p, int K, double PiLow, double Ubound, double emax, double *pars, int lim, int resN, int sph, double *Pi, double **Mu, double ***S, double **OmegaMap, double (*BarOmega), double (*MaxOmega), int *rcMax, int (*fail)); void OmegaBarOmegaMax(int p, int K, double PiLow, double Ubound, double emax, double *pars, int lim, int resN, int sph, double *Pi, double **Mu, double ***S, double **OmegaMap, double (*BarOmega), double (*MaxOmega), int *rcMax, int (*fail)); void dsyev_(char *JOBZp, char *UPLOp,int *Np, double *A, int *LDAp, double *Wp, double *WORK, int *LWORK, int *INFOp); int vecMin(double *x, int p, double (*min)); int vecMax(double *x, int p, double (*max)); void Anull(double **X, int ax, int bx); void anulli(int *x, int p); void anull(double *X, int p); void XAXt(double **X, int p, double **A, double **Res); void XAXt2(double **X, int p, double **A, double ***Res, int k); void cpy(double **a, int nrows, int ncols, double **b); void cpy1(double ***a, int k, int nrows, int ncols, double **b); void cpy2(double **a, int nrows, int ncols, double ***b, int k); void tA(double **A, int a, int b, double **Res); void matxvec(double **a, int arows, int acols, double *x, int xrows, double *y); void multiply(double **a, int arows, int acols, double **b, int brows, int bcols, double **c); void cephes_symmeigens_down(int p, double *eval, double **A, double (*determinant)); void cxS(int p, int K, double ***S, double c); double qfc(double* lb1, double* nc1, int* n1, int *r1in, double *sigmain, double *c1in, int *lim1in, double *accin, double* trace, int* ifault); void genPi(int K, double PiLow, double *Pi); void genMu(int p, int K, double **Mu, double Ubound); void genSigmaEcc(int p, int K, double emax, double ***S); void genSphSigma(int p, int K, double ***S); void fprintOverlap(int K, double **OmegaMap, double BarOmega, double MaxOmega, int *rcMax, int fileN, char *overmap, char *overbarmax); void printOverlap(int K, double **OmegaMap, double BarOmega, double MaxOmega, int *rcMax); void fprintParameters(int p, int K, double *Pi, double **Mu, double ***S, char *PIfname, char *MUfname, char *Sfname, int fileN); void printParameters(int p, int K, double *Pi, double **Mu, double ***S); void freadParameters(int p, int K, double *Pi, double **Mu, double ***S, char *PIfname, char *MUfname, char *Sfname, int fileN); void fprintData(int n, int p, int K, double **x, int *Nk, int fileN, char *dataX, char *Nksizes); void array1to2(int a, int b, double *y, double **x); void array1to3(int a, int b, int c, double *y, double ***x); void array2to1(int a, int b, double *y, double **x); void array3to1(int a, int b, int c, double *y, double ***x); void array3to1(int a, int b, int c, double *y, double ***x); void setseed(unsigned int s); long genseed(void); double runir(double a,double b); double rnor(double mu,double sd); double rgamma(double alpha); void rmulti(int n, int K, double *pi, int *Nk); void rMVN(int n, int p, double *MU, double **S, double **X); void genData(int p, int K, double *Pi, double **Mu, double ***S, int n, double **Y, int *id); #endif /* OVERLAP_H */ CARP-1.0/Overlap.c0000664000076400007640000004711111357671640013615 0ustar volmelnvolmeln /* CARP: Clustering Algorithms' Referee Package Copyright (C) 2010 Volodymyr Melnykov and Ranjan Maitra This program is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version. This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details. You should have received a copy of the GNU General Public License along with this program. If not, see . Authors' contact information: Volodymyr Melnykov Ranjan Maitra Volodymyr.Melnykov@ndsu.edu maitra@iastate.edu Department of Statistics Department of Statistics North Dakota State University Iowa State University Fargo, ND 58102 Ames, IA 50011 */ #include #include #include #include "array.h" #define Inf 1e+140 #include "overlap.h" /* computes parameters needed for computing overlap * p - dimensionality * K - number of components * Pi - mixing proportions * Mu - mean vectors * S - covariance matrices * li, di, const1 - parameters needed for computing overlap (see theory of method) */ void ComputePars(int p, int K, double *Pi, double **Mu, double ***S, double ***li, double ***di, double **const1){ int i, j, k; double dtmt; double *m1, *m2, *Eig, *detS; double **Ga, **L, **Ga2, **L2, **Si; double ***Sinv, ***Sh; MAKE_VECTOR(detS, K); MAKE_3ARRAY(Sinv, K, p, p); MAKE_3ARRAY(Sh, K, p, p); MAKE_VECTOR(m1, p); MAKE_VECTOR(m2, p); MAKE_VECTOR(Eig, p); MAKE_MATRIX(Ga, p, p); MAKE_MATRIX(L, p, p); MAKE_MATRIX(Ga2, p, p); MAKE_MATRIX(L2, p, p); MAKE_MATRIX(Si, p, p); for (k=0; k (*MaxOmega)){ (*MaxOmega) = OmegaOverlap; rcMax[0] = i; rcMax[1] = j; } if (j < (K - 1)){ j = j + 1; } else { i = i + 1; j = i + 1; } } } if (asympt == 1){ while (i < (K-1)){ if (fix[i] == 1){ if (fix[j] == 1){ s = 0; for (k=0; k (*MaxOmega)){ (*MaxOmega) = OmegaOverlap; rcMax[0] = i; rcMax[1] = j; } if (j < (K - 1)){ j = j + 1; } else { i = i + 1; j = i + 1; } } } (*BarOmega) = TotalOmega / (K * (K - 1) / 2.0); for (k=0; k eps){ (*c) = (lower + upper) / 2.0; asympt = 0; GetOmegaMap((*c), p, K, li, di, const1, fix, pars, lim, asympt, OmegaMap, &(*BarOmega), &(*MaxOmega), rcMax); if (method == 0){ if ((*BarOmega) < Omega){ /* clusters are too far */ lower = (*c); } else { upper = (*c); } diff = (*BarOmega) - Omega; } else { if ((*MaxOmega) < Omega){ /* clusters are too far */ lower = (*c); } else { upper = (*c); } diff = (*MaxOmega) - Omega; } sch = sch + 1; if (sch == stopIter){ (*c) = 0.0; /* printf("Error: required overlap was not reached in %i iterations...\n", stopIter); */ break; } } } /* run the procedure when average or maximum overlap is specified * * Omega - overlap value * method - average or maximum overlap * p - dimensionality * K - number of components * PiLow - smallest mixing proportion allowed * Ubound - upper bound for uniform hypercube at which mean vectors at simulated * emax - maximum eccentricity * pars, lim - parameters for qfc function * resN - number of resamplings allowed * sph - sperical covariance matrices * Pi - mixing proportions * Mu - mean vectors * S - covariance matrices * OmegaMap - map of misclassification probabilities * BarOmega - average overlap * MaxOmega - maximum overlap * rcMax - contains the pair of components producing the highest overlap * fail - flag indicating if the process failed */ void OmegaClust(double Omega, int method, int p, int K, double PiLow, double Ubound, double emax, double *pars, int lim, int resN, int sph, double *Pi, double **Mu, double ***S, double **OmegaMap, double (*BarOmega), double (*MaxOmega), int *rcMax, int (*fail)){ int asympt, sch; double c, diff, lower, upper, eps, acc, Balpha, Malpha; int *fix; double **const1; double ***li, ***di; MAKE_3ARRAY(li, K, K, p); MAKE_3ARRAY(di, K, K, p); MAKE_MATRIX(const1, K, K); MAKE_VECTOR(fix, K); anulli(fix, K); eps = pars[0]; acc = pars[1]; Balpha = (*BarOmega); Malpha = (*MaxOmega); sch = 0; do{ (*fail) = 0; /* generate parameters */ /* printf("Simulating dataset...\n"); */ genPi(K, PiLow, Pi); genMu(p, K, Mu, Ubound); if (sph == 0){ genSigmaEcc(p, K, emax, S); } else { genSphSigma(p, K, S); } /* prepare parameters */ ComputePars(p, K, Pi, Mu, S, li, di, const1); /* check if desired overlap is reachable */ asympt = 1; c = 0.0; GetOmegaMap(c, p, K, li, di, const1, fix, pars, lim, asympt, OmegaMap, &Balpha, &Malpha, rcMax); if (method == 0){ diff = Balpha - Omega; } else { diff = Malpha - Omega; } if (diff < -eps){ /* not reachable */ /* printf("Warning: the desired overlap cannot be reached...\n"); */ (*fail) = 1; } else { lower = 0.0; upper = pow(2.0, 2); do{ FindC(lower, upper, Omega, method, p, K, li, di, const1, fix, pars, lim, &c, OmegaMap, &Balpha, &Malpha, rcMax); lower = upper; upper = upper * upper; if (upper > 1000000){ /* printf("Warning: the desired overlap cannot be reached...\n"); printf("Simulating another dataset...\n"); */ (*fail) = 1; break; } /* (upper > 1000000) prevents nonstopping loops */ } while(c == 0); } if ((*fail) == 0){ /* correct covariances by multiplier C */ cxS(p, K, S, c); break; } sch = sch + 1; if (sch == resN){ printf("Error: the desired overlap has not been reached in %i simulations...\n", resN); printf("Increase the number of simulations allowed (option -r) or change the value of overlap...\n"); (*fail) = 1; break; } } while ((*fail) == 1); (*BarOmega) = Balpha; (*MaxOmega) = Malpha; FREE_3ARRAY(li); FREE_3ARRAY(di); FREE_MATRIX(const1); FREE_VECTOR(fix); } /* run the procedure when average and maximum overlaps are both specified * * p - dimensionality * K - number of components * PiLow - smallest mixing proportion allowed * Ubound - upper bound for uniform hypercube at which mean vectors at simulated * emax - maximum eccentricity * pars, lim - parameters for qfc function * resN - number of resamplings allowed * sph - sperical covariance matrices * Pi - mixing proportions * Mu - mean vectors * S - covariance matrices * OmegaMap - map of misclassification probabilities * BarOmega - average overlap * MaxOmega - maximum overlap * rcMax - contains the pair of components producing the highest overlap * fail - flag indicating if the process failed */ void OmegaBarOmegaMax(int p, int K, double PiLow, double Ubound, double emax, double *pars, int lim, int resN, int sph, double *Pi, double **Mu, double ***S, double **OmegaMap, double (*BarOmega), double (*MaxOmega), int *rcMax, int (*fail)){ int i, j, k, asympt, sch, rowN, colN, method; double c, diff, lower = 0, upper = 0, eps, acc, Balpha, Malpha; int *fix, *fix2; double **const1, **const12, **OmegaMap2; double ***li, ***di, ***li2, ***di2; MAKE_3ARRAY(li, K, K, p); MAKE_3ARRAY(di, K, K, p); MAKE_MATRIX(const1, K, K); MAKE_3ARRAY(li2, 2, 2, p); MAKE_3ARRAY(di2, 2, 2, p); MAKE_MATRIX(const12, 2, 2); MAKE_MATRIX(OmegaMap2, 2, 2); MAKE_VECTOR(fix, K); MAKE_VECTOR(fix2, 2); anulli(fix, K); anulli(fix2, 2); eps = pars[0]; acc = pars[1]; Balpha = (*BarOmega); Malpha = (*MaxOmega); (*fail) = 1; if ((Malpha < Balpha) | (Malpha > Balpha * K * (K - 1) / 2.0)){ /* wrong parameters*/ printf("Error: incorrect values of average and maximum overlaps...\n"); printf("Both conditions should hold:\n1. MaxOverlap > AverOverlap\n2. MaxOverlap < AverOverlap * K (K - 1) / 2\n"); } else { sch = 0; do{ /* generate parameters */ /* printf("Simulating dataset...\n"); */ genPi(K, PiLow, Pi); genMu(p, K, Mu, Ubound); if (sph == 0){ genSigmaEcc(p, K, emax, S); } else { genSphSigma(p, K, S); } /* prepare parameters */ ComputePars(p, K, Pi, Mu, S, li, di, const1); /* check if maximum overlap is reachable */ asympt = 1; c = 0.0; GetOmegaMap(c, p, K, li, di, const1, fix, pars, lim, asympt, OmegaMap, &Balpha, &Malpha, rcMax); diff = Malpha - (*MaxOmega); if (diff >= -eps){ /* reachable */ lower = 0.0; upper = pow(2.0, 10); do{ /* find C for two currently largest clusters */ rowN = rcMax[0]; colN = rcMax[1]; for (i=0; i<2; i++){ for (j=0; j<2; j++){ for (k=0; k. Authors' contact information: Volodymyr Melnykov Ranjan Maitra Volodymyr.Melnykov@ndsu.edu maitra@iastate.edu Department of Statistics Department of Statistics North Dakota State University Iowa State University Fargo, ND 58102 Ames, IA 50011 */ #include #include #include #include #include "array.h" #define inf 1e+40; /* This procedure computes Adjusted Rand index * * Parameters: * N - number of observations * TRUK - true number of clusters * PREDK - estimated number of clusters * trcl - true classification vector * prcl - estimated classification vector * Rand - Rand index * adjRand - Adjusted Rand index * Eindex - E index */ void RRand(int N, int TRUK, int PREDK,int *trcl, int *prcl, double *Rand,double *adjRand,double *Eindex) { int i,j,n[TRUK][PREDK]; double sumtr[TRUK],sumpr[PREDK],sumprsq,sumtrsq,sumsq,discordant, sumtrprsq; double term1, term2, term3; double nij2sum, nidot2sum, ndotj2sum, Wallace; for (i=0;i> stream xڝTMS@ WLvkoO$t` MB/MLpGj;0u 0C{Ző$qf}wo@0Ak➅2dZ!lb7w gYJTxH<_(><ͶuTiGٺ!?><cN\(Ey E6OIUeaA5nֺm!1hB/DZ1YA$Y|0F "h\MR@FԏzM\5yR46oyoy7iݤ˺ŢBuN˪yį"D[=:D;'$Zmtl_틵ԉ>RMG"j 4'Kq'Lz=O'ŪBWИ`(Xdh 1qET V1 #a!es %vr\#(wQA?&"$ӽdyIm)E4fBٮMfxܖ\ ^e>lE/i"9 #:V&$;c{P2PF $]3:j_{ @J}fC;lvkgߦS]eKW|Xu.qŧ3z-}f]]7RkZK2> endobj 1 0 obj << /Font << /F42 4 0 R /F43 5 0 R /F44 6 0 R /F14 7 0 R /F13 8 0 R >> /ProcSet [ /PDF /Text ] >> endobj 12 0 obj << /Length 779 /Filter /FlateDecode >> stream xWKs0Whl-fHO-QSbIi,ٍIb^>}jUG Go@dqT()+X29O֊֚SЅy],j4aC@z"Leoef9FGkm>0NrH BP~,GnqݭEW")pXrCʀ.0Fa<hm0K31N](XFw17& Όc@JT@HN)@I$4嚯*t>6ٟ9cTT\VhSq:TR*i:UH$}ye"ceRjc*bNd*1w|+֓ӛ 6orh IWzm̓opKizW94vl$*/ƻS-zgzX |7aFG.VLG:QEft @D+3& nsi.}H1V Ż)C5g֭'fJwq6޸F,RpPywx]b/@jB};_YAvshw77Ghy)) VzD, QG_ HÄxīxѡE>?˻orN,AIف].7?)25+j--{16 SӁ`B$DŽhtqVV endstream endobj 11 0 obj << /Type /Page /Contents 12 0 R /Resources 10 0 R /MediaBox [0 0 612 792] /Parent 9 0 R >> endobj 10 0 obj << /Font << /F43 5 0 R /F44 6 0 R /F42 4 0 R >> /ProcSet [ /PDF /Text ] >> endobj 15 0 obj << /Length 2207 /Filter /FlateDecode >> stream xڭXK7ϯTeqH7IlS^{R{s)qόRm;&A?|ɃcW7pwss"eRQpwR$y, FlwsQEoѴ(q,qv2aIw'Xn^:U|TG9i0vo߬iKy@`M7"ftMI Gjm[Qnp+c ?h#wZ5s=bSMèYVUvvEaEQBf<B.]\dR#iZlHO߫AWe?1zvͺ^irҵy-[fUm_ G37]N"\1`_mߗNۂzk͠&VEiHqBw#)]MNsnCJ^4Y4z7tVM 7CQ)/*)Mrojo*3B2cm}HIyFÒ sG5UgTC4`=Į^|p! 4v=P܊t G6XB;]VH=Lv;SC™KRƇ FWxaE[p9 #LRUh xk!Y5z :?z/i "`p(t d `q-+I4[R*,P,'c >z6}ƛ~rkKQ7ҁ# ~k; [(D A!Zh{ךb"iuSFmwX^p`V%,疸:rS |2*ID'.nW*—3 Wa h0?FoN}1PCƳZ2}6]q`t\-wxrv)^h; Z9˕]59ƹ-=aurU3Q8O2cI>̷dfpˏ9Ū qzLmsY NϚW 3RPG>1v~%"W"u"W@Dm1%Ԟ9QFwr!R[n,6V]O@>M쮍~tw \,8˳$(AHp[:Ъw7qONKGqju j $gp a ;).ryUa][{@ͯ%x$a'z0E` [z\z>0 Y\ݼ69抩 7f&YrF^ %1d>e> endobj 13 0 obj << /Font << /F43 5 0 R /F44 6 0 R /F42 4 0 R /F14 7 0 R >> /ProcSet [ /PDF /Text ] >> endobj 18 0 obj << /Length 1765 /Filter /FlateDecode >> stream xڽYI6W"6*E]L@-6;Z\-q~}EJ4q^fd}FvgWg/gF \Ζ c;3/Pr=k fދw=&U XXp';5&tv H Yc|2G- ӅE]b8F\Vg1 fFBGiD0ebvr`G5Ak5fX ˅e+6-(6f[YktqbNi-6ȜgJȰzw-32.Ih2^kk$iRÌ#A&ٯ$S[ Iˇ}R&߄_EOhQ܊ZԦR FV=2bx"~&rrԪ`<ŸL$&/&xaTZ/[s҈9gڒxMI004Aɺ䠯p zOu'qg( (dRZx£× 95AN)9Lo+Gd-tH5p'T'~\@b6cRܚTy7jIJ9wG#n]{&yBBC"+îI'ibVCiCff0&.-ÑBT^Ujf-\JDoYC=d̻@dS&EF B4jBY BAm$s 6ԷR W_!1%)Lu]5"ID|a܀q#qa_⚨|cOiTu9uv8_Aa̐GI4>Huj|O^OG2?hmAnc~ #>x9hivYu"n؛Ж6A#Eֻ±~3̧.$M1t7;snzq0Uu]q9j8ˆy |@:ۼ8>\5\I+Av⫡6u?RfغL]@ӄYʞ^_nH?ￃ[`>텎ӄO]{5mz=F`NR`9m$d|Svd.&:GW qE[2!g^Y}+*n-UAgZ!ρ-,yBuTJ/ۥY?Y.W0o{#x(d^0໭oc4Gi=軗NHQ;HsJߝ*jj^A]us*@`iŜΘgne endstream endobj 17 0 obj << /Type /Page /Contents 18 0 R /Resources 16 0 R /MediaBox [0 0 612 792] /Parent 9 0 R >> endobj 16 0 obj << /Font << /F14 7 0 R /F42 4 0 R /F44 6 0 R /F43 5 0 R /F80 19 0 R /F15 20 0 R /F54 21 0 R /F81 22 0 R >> /ProcSet [ /PDF /Text ] >> endobj 25 0 obj << /Length 3429 /Filter /FlateDecode >> stream xڭ[[s~ϯLJD,q}:ͥsIIif>@l1H I.x)vbXp.&g ?qaMZqvDl>LV7ToZ*;jھlj"[o Sl]߮ȣ}QW]oENU$|lʫo%tuJV'0G"\4|s1 c_dM QZGwB2ˮ+ו9juT-wLo&z,뇿, &ELV%` &}1p?%{Ud\WKwY;tjD9XNv&.QU˰+{Yp*\ʺN%s*3x.?]e]su,1vYI)89D'fs=.6\IDky۔?'tyS޵8i;0 ׎wHaA`rlO [pYy>")gU$QYH `u0' w%v+7ln-~['2~ik2j1wh2 6TH5 ~uo˚{ڶ/Q@Xx A[Y6} mNIoeJ,iX b! 2AR[j vTepڣMeglGs[ƈiԮcfЁ m;wx\ fxo SKs`}F ?l[@1zpڅM')~=Vި†$;:)MR+@ʂxOT_9\fn[GdPtK|x"{nIF9j&ؚ@%g!fڥɌI>{W5qNҽaivbr|O5uYFj>lkvz׺l e]F1Ĝ"ZB_+^a/lW#0ʼn?-{D̅>"\o%f0{![`cݗ hk.+`C*J&0~1єۅX ¶pbB2ىX|'GG#z>:A7^F-'7MuT5ZF:-g'7~! (|p0Ėj7ʳ+Kl$ T`ߴ37QWk[iRijޑ^(^Dr{XpA]b񄙃"NA P73KCh1+ / u]>,M3S5•eB-LO7лu1H&_} +aI"d5֟قiۄq,TKE^T\F֛gy6ƍ+$G ƙxE19L?Q8!JXBOaIDc5]젔^1C%0jd&~EB%ٝ ꦧ. L D n;?EPgՃ T.cQjPSTtߏ!Ws.?qWcO)} ؀*e撀@ ޡa7YEDV RuDϦVvTJBUOZ_\0Ez!M<%nw ,$M_=QDYe_(IE*NIlRj4Cx}T=R@)/HÝg 9c:8cܣѣ,؂fG=l) _7rز`a5({/+*D(~L °fg> bȬ{,E;XXIa@w˴X Iy(8*~jnMgD<ʰ  ;/Ǚ^8$ktK&u=4@ĨX f5,k)բPL]pDNZ`)>j8^yg=yCLG4aXwnQa𺮞F2 hq8s]$P!Eu%S.%b1ubsŏq*^.\;btKn¾g=r{xd~qq J@c¯_ ٙg-&xOe"+<2{_&.qab^! dAϔ`uxvShi="ғwg=֛0p%K A`N|1RE"ٓ ކ@Tdxvh .@[FBHs`A|}&k endstream endobj 24 0 obj << /Type /Page /Contents 25 0 R /Resources 23 0 R /MediaBox [0 0 612 792] /Parent 9 0 R >> endobj 23 0 obj << /Font << /F44 6 0 R /F42 4 0 R /F15 20 0 R /F54 21 0 R /F43 5 0 R /F81 22 0 R >> /ProcSet [ /PDF /Text ] >> endobj 28 0 obj << /Length 1692 /Filter /FlateDecode >> stream xXKoFW>͍m HrĵVUߙ!E*tbrwf|v7_^pzx3!Xj-ngEiDcPչTiVӳ^B][Sڢ?/~kZ!%ʀTwkȉh$[LrѬf,wI)9G:+'xv7E9Cęf|(6zN;7mFx,TF4JgL_T cV6ֿhf頜hԁBYcϫW$qsLHL $g0Ipy?<6{i{^L`Wۚ8]?_55Э͋fg>SmYtޯ}=I`ozǫp vDgm6ȦQN dź@ecع0Ah;L!f%04Ex%umVyQ-Z|:\Uڂ~,^x;tГ վKdumW(PZrsnp6T}x jΆUsisa7'L^fH& ȚKNuvg)670C%< b۸x5s럐ll}yƒ.nZ$-~Ӵhut˟-..gm_]]nX3İM;`[|B&,I 5!SEBB57@k?>PB 2> 6"zLDIvCI;&fAgHqǎ䱜QdzXI+x:#ƭ@77w!$SNs8 1 Ds9hGh Gx!Bt[lf~8VUN«WK*[[w0+{q-1 U$b?Kq%L]UV49`:$w6;PQOoT46a5>x uPAB[J. [.6t;nɚ5U&8@bXql@6@iEQ?xcb%8PtPyվHZB4.Nߥ%G&̃PTa6T,孻6ۋ'q6;IxۦpnŗꞪ_ .n-cxqHƇX@e,6E[؆zi$8ӇlF iwɛ-~(C|->%#F! 2!v(EdP*f PnJw$,A2V1'x A*H 4A2I*UȐߚ>䩂 IiqDAAI*I똫 |p=w:f4b{[!S&D˙RְĨyRW&qA*tL='܁!Ds_P@Fh'3G~z#3]lݜM /oDԒ1ux pmMшr_n endstream endobj 27 0 obj << /Type /Page /Contents 28 0 R /Resources 26 0 R /MediaBox [0 0 612 792] /Parent 9 0 R >> endobj 26 0 obj << /Font << /F44 6 0 R /F42 4 0 R /F43 5 0 R /F81 22 0 R /F15 20 0 R /F54 21 0 R >> /ProcSet [ /PDF /Text ] >> endobj 31 0 obj << /Length 877 /Filter /FlateDecode >> stream xYO0篈P;MiiID-j6S6-i (Iũz}r?_y2Lxn`k@b,&80T(k5"'ǽÛ,O K>jۭ~4 p[g[ n%n?@)ZҤ* 3֬xBP^X% nBDBj]b(b|8J7ݩ*B+h:1Ck,'S0zKgQXkt 2<_0z0:5_> qq .^ZEVnyyήn2/E7+㧟~}4aڰ_/Ӱ4iiZf猱L}|?W$>YOse+^keeM}=&A4/ήҢ/'4=ͯ2ϊa`KGO<$.S;ǟBwrj"?% iES==UZ?97ikC-ml߮ endstream endobj 30 0 obj << /Type /Page /Contents 31 0 R /Resources 29 0 R /MediaBox [0 0 612 792] /Parent 32 0 R >> endobj 29 0 obj << /Font << /F42 4 0 R /F15 20 0 R /F54 21 0 R /F44 6 0 R >> /ProcSet [ /PDF /Text ] >> endobj 35 0 obj << /Length 826 /Filter /FlateDecode >> stream xYo0_i/66uژ4!FKIئ3؅ > t8>~wAADwGӷKB`ԹA)"5ъE~t}LO9Qcgz#hR'-Uy԰ρϜ-j I$WQ 8\Tz!# /lTN(+I(8= q`8&T 6eZΓ"MqҒN; gya2M~:(.nMOƈ}el1/>+ TDit;~QSU` $$X   cH[6Rf8l(p 5Hk8 5U80BM&f2镻Bkj r!6 Lćt4ѭ˽kSo2Nf\hNd[l%?,L,ڷ.HEu5P$fl+b>^$m<#BYt`66R_M)A`J2[s7[ X\BW۰` X12=9W2ͦ:QkVcaq4>v@pXxa]i|O}ǟAMZ£Ot=1?GV MmOi];uUeE endstream endobj 34 0 obj << /Type /Page /Contents 35 0 R /Resources 33 0 R /MediaBox [0 0 612 792] /Parent 32 0 R >> endobj 33 0 obj << /Font << /F15 20 0 R /F54 21 0 R /F42 4 0 R /F44 6 0 R >> /ProcSet [ /PDF /Text ] >> endobj 38 0 obj << /Length 890 /Filter /FlateDecode >> stream xWO0篈TZ?-PHm&m{Zfj*q?v섦2(wvg߹wZZzwY|rzB (@膕E5'HդyV;GK :#;n#d#=Jtr㊧YBjRU/@P!!?v,ERQ/l!y܈Ww7эGw{²\Aq^XӽZcJ`c@ |;H́A= endstream endobj 37 0 obj << /Type /Page /Contents 38 0 R /Resources 36 0 R /MediaBox [0 0 612 792] /Parent 32 0 R >> endobj 36 0 obj << /Font << /F44 6 0 R /F81 22 0 R /F42 4 0 R /F43 5 0 R >> /ProcSet [ /PDF /Text ] >> endobj 41 0 obj << /Length 1551 /Filter /FlateDecode >> stream xXKs6WḦ́MδI&-́!-*In}2"MqJp.~yA͋W_ \Dl٭[EIB6K [@) ަЯ!qu1j, űt3b_7ywUn}MZEݴy]u5mua*YH½EݱYCh$*2G+4 rL2a 9܈,kl@&hX(u,RW-k.u82Eo[ Ȭ1zd6yZm Ÿ"L&Zi$P3 ԇ.U힅n*1J n&tuϩ#qN)Y/}mׅ۪wޮW<9%2Țl8B@퓃muĮKG9sI9k,.}6XUi5Ftߤ/,vZX\庱0DX֗6UsKg iϫ̕-YMefr|WEms3FM%y=OyE^qFaI ߄'* KcKʒ`YL}yv_a㔣j'.|33ʢ>I\İ x(6'1ջN%bcHz|]X(8TϹпC%!y=%p>_l8"`\v55C4CMC wĆL PNuF#T4nQ'1[ ya:{dn^ޮ?ǮCV71o& EK.2рM|kvwj?ۨRUsl &(H{؅s. F3W^86ܸ1˺s,o{95$xhƟ/6, 1?GȘ T "8P>;"W4;4sS%D, D4F/(as5f~^%,o6Mk'% FAg- #tή?Y ⨉<0;98#1N2ES?yy/! endstream endobj 40 0 obj << /Type /Page /Contents 41 0 R /Resources 39 0 R /MediaBox [0 0 612 792] /Parent 32 0 R >> endobj 39 0 obj << /Font << /F43 5 0 R /F44 6 0 R /F42 4 0 R /F81 22 0 R /F15 20 0 R /F54 21 0 R /F14 7 0 R >> /ProcSet [ /PDF /Text ] >> endobj 44 0 obj << /Length 1268 /Filter /FlateDecode >> stream x]o621$EZb=3kvAX-;(4 ҭx;À~~G~cDIQFHFJDy6 a$~!D~y~-!I@ pc$Jt]3;no$O$ cc"R;eO!5IDkv}{#֛^Jd{zz *mZ7pDi7+] 7uןpFW&>usZ=BWi3B+\ 3ajF1a4FqyaG)ȋ04@MLF$DJW+*וJ[_aڈ8b5_`ʠjVmU>TiD ` srNG-*_]-j5,rcIKo:}֧-kK)_ Z>RS\EL9 ^ }]3m30'Az_ ?Z6 if:#:'ǽś2ɛ_W|͗D*|;e6o t|Lf1ZH* QI ! &Of\ &߬a׻)''v, Jq Fv|[c5[VkhV 0#DE.b0r]ٺaːs"my"]ۓ Meۺrl~LFxp0E$9IWG4nr endstream endobj 43 0 obj << /Type /Page /Contents 44 0 R /Resources 42 0 R /MediaBox [0 0 612 792] /Parent 32 0 R >> endobj 42 0 obj << /Font << /F15 20 0 R /F54 21 0 R /F14 7 0 R /F42 4 0 R /F44 6 0 R /F81 22 0 R /F43 5 0 R >> /ProcSet [ /PDF /Text ] >> endobj 47 0 obj << /Length 1043 /Filter /FlateDecode >> stream xX[o0~cf\ƺv.^>bIP.?v줁u(V>s|0rrޟmtv~͘1 9'Np29^0N8\ Gчnu[0RY&ޡd>QoNek yB.l'!kdEFT#@ryYbtHd*RY6S ]*3mUvâSӫZ%X0rzm7|GUȋN hƲ&:|&/*eJfT]qCYnD-.K:(R>+E0"܍1zscm^0=^lhƋ!f[q` brx B`1RզCX+b8`10I@#R"1Axc4饢qݕ/,Tu\r!MʹtZVR\̇biC(Tl;`2bW2Du,5Fց 3t:u?<AGHIT$D@x}IOHG' $j$!s~ H"'&'$ᣓDIHC=I$$|\\խ/}zgA-c3d9x`MoeZT| @vT5DMѧ;T^][CJUW7Kjy@uYe^j*{7<_IQ4PZjPԙ6e8ojeڧ0|)h܀zχuYJWV(T)ZZF\ra'`Å.!*r{6^3a4Noo.bGlOM[[L;>UCf[yD&khyrw]uyryIdA2r ?( x0w?;եQ;P =2o9'qBZOu3Bdb|OZy9&cʼ×if?}N>_(suqx[&;[_2 endstream endobj 46 0 obj << /Type /Page /Contents 47 0 R /Resources 45 0 R /MediaBox [0 0 612 792] /Parent 32 0 R >> endobj 45 0 obj << /Font << /F44 6 0 R /F81 22 0 R /F42 4 0 R /F43 5 0 R /F15 20 0 R /F54 21 0 R >> /ProcSet [ /PDF /Text ] >> endobj 50 0 obj << /Length 2156 /Filter /FlateDecode >> stream xYYs~ׯb^U8(؊U` _`EF{{>76V3rdr+rW۝Gؘ7۝-N/9L?ADyivR8KܔD]x4J4~HlȎ .434]BTIJ,+p[_5MHTQYAe~_xlT@ %BF…wyh& U4V&deew8zC3\;y] ;D|>^|Y .>O??MOerƉ]TzTZڥpv&㫻G.A.Tt\WhWU+Թ@IiE֍k"<8,"p& _d:Vbf л0$>q9_nq z콈;kbl= @b|lkJx}wq$#S< /ґ9#f؅A5԰%gtCkϣ>ywS7lh{&#Vݙ?Ԧ ou'>LR_6f O9 D)LJ~XI>aj`9 NdlP[t%WuI^˺X^$kqIg~Y'.K 8c|b2p_w.vF.}gϞ!(C‹3G俹#\wT%?@w@"ϐ/UWY#X[VO 1>IXWƪXUI >{yV2ՈD2<+$0+g !]>{B:*Xbbt%utf4Xz\1ω]T˭tWQr&6RĐ-R)eU`F_MJ)JBJ/LޭL-*$𺡖j=]5k:J":ZRGc@¹)O|lyG@f=uj=>b 'L\3/IeyK$)?~LHͱ͚y95 jb#L*,?|^;Zl^aΔaelv_!K.9,dT2F>wCT|wec-$Re3$ !A~bK endstream endobj 49 0 obj << /Type /Page /Contents 50 0 R /Resources 48 0 R /MediaBox [0 0 612 792] /Parent 51 0 R >> endobj 48 0 obj << /Font << /F43 5 0 R /F44 6 0 R /F81 22 0 R /F42 4 0 R /F15 20 0 R /F54 21 0 R >> /ProcSet [ /PDF /Text ] >> endobj 54 0 obj << /Length 840 /Filter /FlateDecode >> stream xXKO1+VbQJJ`J׈8) u9:jS,=X0nO?/~m&̒G > endobj 52 0 obj << /Font << /F42 4 0 R /F44 6 0 R /F43 5 0 R /F15 20 0 R /F54 21 0 R /F14 7 0 R >> /ProcSet [ /PDF /Text ] >> endobj 57 0 obj << /Length 1151 /Filter /FlateDecode >> stream xڵW[o6~ϯУT I`Yf+֡H}-V#K.M_CJـ0˹;߹H.ϫL%VۀQFx&L)d6e_2K{ZrIM3@pF$30M WUT樖S8D*,%irwЋ^\Nуs"tAH;rJYkT&i"#Y*hT^}RyoTxGl #TQ ľ A%DqnuK$BXo׻Qi>QQޗwjz+1n'\5[Sngj* XN88-*ܘ L1)hUS[$eSw6}6^xF$Sh R7wzgHHͧ% #.;tqm1=}T(+LE=4l f\Gu5辬wiM{t٢ʐe:Ў0U=Is5O,z0r5IP2~H'Jy@2ԅiI7}]=tTQJ}cfCו]1O"AY[ͦ <LjDi % 'DV@za8xϺlΠ6ްG`O,eX \ZFC ma'6M;EA~5.[.NQ{>ALud/*q|5 m+׫ RhL?π Lޡ%3鯘,RI* )P ?oݐS#2+b\0wT}dŠÇa&^.®G/\ r6Г8&c1;gٰ[Hd9bk-!lú*}EV . nޙ"6׭%R‰>}}Оŗa& _ESŸW f jt^3|03 M_1uu/CN8f endstream endobj 56 0 obj << /Type /Page /Contents 57 0 R /Resources 55 0 R /MediaBox [0 0 612 792] /Parent 51 0 R >> endobj 55 0 obj << /Font << /F15 20 0 R /F54 21 0 R /F42 4 0 R /F14 7 0 R /F43 5 0 R /F44 6 0 R /F81 22 0 R >> /ProcSet [ /PDF /Text ] >> endobj 59 0 obj [500 500 167 333 556 278 333 333 0 333 675 0 556 389 333 278 0 0 0 0 0 0 0 0 0 0 0 0 333 214 250 333 420 500 500 833 778 333 333 333 500 675 250 333 250 278 500 500 500 500 500 500 500 500 500 500 333 333 675 675 675 500 920 611 611 667 722 611 611 722 722 333 444 667 556 833 667 722 611 722 611 500 556 722 611 833 611 556 556 389 278 389 422 500 333 500 500 444 500 444 278 500 500 278 278 444 278 722 500 500 500 500 389 389 278 500 444 667 444] endobj 60 0 obj [272 272 761.6 489.6 761.6 489.6 516.9 734 743.9 700.5 813 724.8 633.8 772.4 811.3 431.9 541.2 833 666.2 947.3 784.1 748.3 631.1 775.5 745.3 602.2 573.9 665 570.8 924.4 812.6 568.1 670.2 380.8 380.8 380.8 979.2 979.2 410.9 514 416.3 421.4 508.8 453.8 482.6 468.9 563.7 334 405.1] endobj 61 0 obj [380.8 380.8 489.6 761.6 272 326.4 272 489.6 489.6 489.6 489.6 489.6 489.6 489.6 489.6 489.6 489.6 489.6 272 272 272 761.6 462.4 462.4 761.6 734 693.4 707.2 747.8 666.2 639 768.3 734 353.2 503 761.2 611.8 897.2 734 761.6 666.2 761.6 720.6 544 707.2 734 734 1006 734 734 598.4 272 489.6 272] endobj 62 0 obj [600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600] endobj 63 0 obj [585.3] endobj 64 0 obj [777.8 277.8 777.8 500 777.8 500 777.8 777.8 777.8 777.8 777.8 777.8 777.8 1000 500 500 777.8 777.8 777.8 777.8 777.8 777.8 777.8 777.8 777.8] endobj 65 0 obj [556 556 167 333 611 278 333 333 0 333 564 0 611 444 333 278 0 0 0 0 0 0 0 0 0 0 0 0 333 180 250 333 408 500 500 833 778 333 333 333 500 564 250 333 250 278 500 500 500 500 500 500 500 500 500 500 278 278 564 564 564 444 921 722 667 667 722 611 556 722 722 333 389 722 611 889 722 722 556 722 667 556 611 722 722 944 722 722 611 333 278 333 469 500 333 444 500 444 500 444 333 500 500 278 278 500 278 778 500 500 500 500 333 389 278 500 500 722 500 500 444 480 200 480 541 0 0 0 333 500 444 1000 500 500 333 1000 556 333 889 0 0 0 0 0 0 444 444 350 500] endobj 66 0 obj [556 556 167 333 667 278 333 333 0 333 570 0 667 444 333 278 0 0 0 0 0 0 0 0 0 0 0 0 333 278 250 333 555 500 500 1000 833 333 333 333 500 570 250 333 250 278 500 500 500 500 500 500 500 500 500 500 333 333 570 570 570 500 930 722 667 722 722 667 611 778 778 389 500 778 667 944 722 778 611 778 722 556 667 722 722 1000 722 722 667 333 278 333 581 500 333 500 556 444 556 444 333 500 556 278 333 556 278 833 556 500 556 556 444 389 333 556 500 722 500] endobj 67 0 obj [600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600 600] endobj 68 0 obj << /Length1 847 /Length2 2667 /Length3 0 /Length 3255 /Filter /FlateDecode >> stream xڭy<{I)6c{I(31Ecdd e$ !K==j,E ![:orsy|?u݊6v&;hN2P0b-H}!t ѨD]]-/@jp4JBASwDdL H'qTcxv< $Ld A@DdI_JA 8 ! 4*9 5f_`ܗLQ~ߝq9 ۗ,ҩ?GςaAɗsՒ#&T"P>ǜlH '#>H%Lž.sV߫-HTƙoՈ4{Ht?vur13*F Qt\nD @@gàT}`O&!?`V?]P OPp;8;NR ]'Ͼ^3𷥅`>dc:~ԑ(66&GKT_ڃx` q1:Yb^eٌˆ =ܨ[: xxY{QbIsGy 6ʦNv5]#PcVt%(|>I/7 Ov~u~0'n~0W7ԊV՜5Ĩ`wJ*zF{:ru/#?n8Rl&Y:/y`ۥp< |rtTa |w@B}sX TO?Ð%yoAsI/,^gEXȢ3ɰf͡)<1奉7ArKDg7АD~![W=f_4 C+,h9:Zܸ:]hv3gN 3ZhȖseVBS+"ZƲT}џRlQ1q1Os˃}pmoUo'Sw[}Jkm+ÑXzo\ =q aU$:~욠fK}#U.U71r ^ {}9 Vjqti$cXJ-@3)2ך#p3ڷpz]\$9ﭏ{ύњDfCݮ rDe M{)PxnD5eK'%4Lc {18vfcɌF1PԤ7{ՎSGՑ'2 Ւ?Uj; pՄOhW^]ɡ͙'rN~gU[W! 56h|d]I +4ۀ˕,h+;ϧGARohslC&]q-N#X*ZL,@2zOV'kUm(X --Jؗݳ2n}ɎsTE6:Ry76*+,qӘYxӌ:YD9_ݟiPj74Kz&&_f`DpƼ$Xayjʠ ?./],\LniR?w&sݸb i x⾰$h {`;t4GY8x̍C`VWU-z$ړK"Im+|7C]'}Մ5*0ia[[bW6Op8fV֮"VᲴ,_5y{ۋs"^ۼf,ɶYE&w:=aȤr4SCG f;4uxazdTй˜&7C,0y11|IoY"gߪTu{=m%q: s×<>Ã>DR GzE]=KάEǻPXZسx͕ =3$FWgŇΞG>ؐ,( qɨڪ_ ҵx{CqeCS6H=h]r\pgogKk$1 ;&7Pgtun16vjRv>uҴdߴ8D6-f6▗"dyNWsJV ~3N[unZ|]팘FNm>J*XR oӛ׍Q2-xO$aҤ Es"O]]-,0^No?4 endstream endobj 69 0 obj << /Type /FontDescriptor /FontName /FGDBER+CMMI12 /Flags 4 /FontBBox [-30 -250 1026 750] /Ascent 694 /CapHeight 683 /Descent -194 /ItalicAngle -14 /StemV 65 /XHeight 431 /CharSet (/K/X/comma/i/j/period/slash) /FontFile 68 0 R >> endobj 70 0 obj << /Length1 1003 /Length2 3085 /Length3 0 /Length 3743 /Filter /FlateDecode >> stream xڭyeT)Oht[<˜2ZMe %ζ-ىeb ͙:z! 2kh$s&񋭎?u洁b]t(~agJ3;FƒO[Ex)Gvֽ4QTU| L_vCJu'Ps ]n4➓}Fǯ/!/C se ?22`jcK h6ZF lȘ̐rE Rn+nyE7䓌't]#MRy11lL139@1o=|ց8=Ǥ)_WF¾ʪcOTsV'*EJZuR9WB RK^N,zx4N/9S[)10G -8${j@Ug{Ld=\Jy$.{/(ٸ^jPvLẫQ2)#>j7FPbÒ]!2ES&-H |t|̑O5cUA}+թ3WnXY ;V.fzC59rWЙwyS8|ڽcK6!Kѻhu(2ra@J[]a8o[#/?:L#J~ڶ)E*OrpTNoqFٜŨM> }u~G[C.] we?5NV kA߮˝Ud.|K-=&ǴޒzRt"iɍɽ2]$- FDt^]Pz"I,ݽiQ»Y~)vki]Kac[ gbdǃQiG4#K}>ih4"oٴRXcr}o meGtZdfwMeˎ_S/) YxZRyc&omr+vvWn;:hBq2()b1 rM9r[\u)1 *3;Qy1r$~O}D"%mzyD) TdM@щؓθm.{˗̶V\TD2oiC OMy`TS@G8<Y i` +E7Nڱ}:i{<$ 5oq6MM6?7\bOmm*C*lFʽ˯ ;*}./<{]&(~Zx5U?TW}XdzkoahuNͨ<"1qӋm%mzM/XY% $Pbg[ l?Yfw4uȗ{hz.51Odp&Q![rVVg9q~$;KQY!5~nkO l`&b0}+ѥ֔cy]>'4qV屚5aᘕ>$Y_Co˓/ţL̔{҈1lbS^9s蠅G vi|Y=G}=jkv#j vw}%@[/s)z՛[ϫK0 h:u+>p* X5qKe|C]<+;g>b(p:>kIڡ:YgjDpF u]ŧNt &+YP* ooI"-'7ЖeM; G3ç ~5=g\* Z"~*f mXͼ33cOWbqLRYgWL9?> endobj 72 0 obj << /Length1 811 /Length2 979 /Length3 0 /Length 1535 /Filter /FlateDecode >> stream xڭRkTLk8) ӽ4LI(# Y3m={BQT*qP1.+K(SBg:{{Y2B!ށ0 dz(LaS(@sy &&ILA;o(IaT{HaR Y@P/T`1Bu(b@@0)$G0JCIm &m!pE4@P-"PajhsZʁ)JLEA(J@AI|t\ `jP֟T(pXWgS 1 a4` Qjo[_\ é%Q(|WN4`@>_oC>ťru0I z b vfE?t4@F6TT"A^3AVbZ pqlF؏)xygN wqc#IũEZ0:ՠRN\uz[Ioa͑lcCr1$Vl-=/h,K2=9@:{Yj"ڱѲНb1.8_I66I9#^kiZ^xhtSx? dYs/o>Xd&ޖ5tuo;?\mk\?5]rj$+bMWrԘz])]2*>Jndm`R3E(ʜ 4ti5!'CI Ž/#dXtlX{ܜXւa2;72WNpPtcTK]fyDj;ۼۍCMÌ-m ?HY G\jijht-c&+~Z)lUAlIȭsoXھ'/y.ҹ1,Gğ#8Sy];gAɨ&]a՘ds Wb-2y2x/32=_C:SN ~j98,Qd[=iDM{gZG]=M֋=d%{=}/WZu@> M>/K/tG?/{Y,8wyzMb31mlQúF~ o Wi^ &T,Fpee9ȍN+]ޚ{gR-cKz =#w3;tcվF& wLq*|_ntsI# 6I2nk2eozIcOiX]OE g kl0+Lo5_ID:Aޝ4q'N>1kԴt3V6Rmջ}C-%o{3`yb7WZcߴU6;W,8^ #˒FH76lgV/hwXS^|Ur,.xz endstream endobj 73 0 obj << /Type /FontDescriptor /FontName /NWIAKC+CMSY10 /Flags 4 /FontBBox [-29 -960 1116 775] /Ascent 750 /CapHeight 683 /Descent -194 /ItalicAngle -14 /StemV 85 /XHeight 431 /CharSet (/asteriskmath/bullet/minus/similar) /FontFile 72 0 R >> endobj 74 0 obj << /Length1 752 /Length2 721 /Length3 0 /Length 1247 /Filter /FlateDecode >> stream xSU uLOJu+53Rp 4W03RUu.JM,sI,IR04Tp,MW04U002225RUp/,L(Qp)2WpM-LNSM,HZRQZZTeh\ǥrg^Z9D8&UZT tБ @'T*qJB7ܭ4'/1d<(0s3s* s JKR|SRЕB曚Y.Y옗khg`l ,vˬHM ,IPHK)N楠;z`9;hCb,WRY`P "0*ʬP6300*B+.׼̼t#S3ĢJ.QF Ն y) @(CV!-  ,IH 3NNպkij`hdj`naT0(5\Ae15"5d떬Vֹ.U牵/oQ7;36t^bzɒW m;,}Q+zp^/lZ,bB~]xi<\6cN#iEgu*\g^ ߻-4IKI bw_ ۟g{ 4$=?ؽlrڋ-СCg~ЖԭVw\S(ȸ͎'q9.qblSiqkx(}'ax9vUMyQ=)zvf;zBO7L"w+X*8fwvܲY&{z!g*_ Xmu0E5ae_پV/W*m˳ ʖ%z4U< ^}ru/z2){ hz l5kmio ~[W.l5əq3>ii7RA,+5zI:/zcfs}Ycow$ǯN-m欙eXŋ7;:FyXCҷOk24oh_rcٝS޺^?A03jWŖʛg^rdX#r9_jr4r,&fzflN (YC endstream endobj 75 0 obj << /Type /FontDescriptor /FontName /SDVCAT+CMSY7 /Flags 4 /FontBBox [-15 -951 1252 782] /Ascent 750 /CapHeight 683 /Descent -194 /ItalicAngle -14 /StemV 93 /XHeight 431 /CharSet (/asteriskmath) /FontFile 74 0 R >> endobj 76 0 obj << /Length1 1612 /Length2 15808 /Length3 0 /Length 16640 /Filter /FlateDecode >> stream xڭeP$ܒ-[. ָ[ָCNָ;37Q3s\+wFTDQ(2%lA ,̼ #g'y;[93௑Bh7y@8 `Aٻ;Z*4tte'`7,lmU@0eR )-ldma 0:ivv&K `p[t3m,~X8 mA{X[;C_FSs9;Z؃*K'On'nH;cJ/_/ eX8[_4,l=hfhb tr WzC{{kݶWr9MXX4mfaϠؚXm7q _ gfh04vM@STf@G'{tV0;^0 gX: m,{&$?pd@!bkWfF-$-܀&J csNˮnktU_0037? hk_ԙĴTNW_AjRL d`epM/:-KffW_'#aklgϬ mMq;;:U_/oqנn@c%;c ˴tP-NЄn_ P}iZQ_]oZOŏ`)ޏ6So#=T)d4h\tLHgQ rP:'UK`.i\ 0}Sc:Pk O(F{ rb)\(#rt+V S;G-bU䩬0~X(  cyO󤾎VcVCȐfxc _#-H%/j }"g:ї Y o|"}cٚJrh?$}Rcx ۈ 4 t("=k@tN:]<^|A^v^UZg"m1~2џpKmw`s@'[* tum}EݓsW2ONE3OVq qܰ٠"dX;9oPiD\堝qSn6އmullv/s Jӎ-(_G, 9N[R7b=JG DCjp.4lS_>G'C,!ľ$N;hba[i!S73=TmNzXٔ4xX8 M rD ORN9 usQ9:[A6wp>-Y[?'9PU*efzO4uB\Iס| +Ղ*u%lo8A&tv*msO9_taaS8#:9=ݔ$Sr/V9z`,s@j疋;Z\s^YJдJN;8}b^?ഀqhIV=U84,gA5Gɜ?竑m-g*13Se5v!bI@yϟ-y oSi4s-IȅsWVξRG% |..dYIgy/<6q.?%iUy}.b ^WzĪ(ڎ`|T"&K@8,OHm;3L:h%J:3.\)n|xQx ]Mtn.k< 5|R2!'ƀxZi$J%HLad}cMbo]#O6TZ\[L;:1f4@iݶ8jmQ_STղ RL'K[HF* lv_Ӕ" 7xo&9 C<t+c~Կy5?lw!B CVpg臾snW@!:U],rz̲sHuP '7VDoeJ2nUy9Ƥp^ʬ ~IoHcoqCScNJ:/Ҿe lB~V !3)f<ʙ8OK9:lwDw)mɊ?0RX{x2꾻]3] n+ O9*d' t(_`xPf_[,bzIucCJ,"boUKӵ$zo?}VDnFGӤ崫#6K'm։u6cd5[k}8,tͽOeS s_E bF_/QO nzH]d.EIw.ԈaRgQ_PK(]zC@t:#/y ñpw{łBIakǺ'A |kW7_!$ d\&IƉ8`U\jӎb3> x)4 Yv(2^TKnBv/o QaE)ў]ʨ- eOE"f9J0l.FmS{优;u. eR㑭 #vѫ! lQZbiRuSa {H6p{K)5i) +ڞVG'-*vI>SgQ!׬L8cyQRDjv`,9ӂ!VQSAs\4:Ү&݄ԿhFzA^Bc)ɠ3kG#LfI2|tҔzl1.*U.7jL1W:뾤Dnb/'XH~%H& k\647f8-֒T\jXë_I+y;(dPVBbMơlYE ;L|Ԙo[t ]KDŠs ~ⶫ77x ~.Y/I;*ۅY\QǗbdLDq8],\B' ?aPT)شA|2O5xZ#Fqܗa-07 o*y˅ Ve(VpX[|X#),πyR ou-7 ]N䇞pRQBzaFhe`$03>ϩ.=#&nGv 3WVBtN#m}?-v]G%urJ M?^R FS{+F OmdX^>c熈=BI5F`Ys#>zo~a a'J_BRUZSVyt KOms"P2yr*YJ7jrj9PF(yUGCfΏSduyy! Շ^`1??TٚSU(t _P/m>C9glpTlHNط#A6v#Ӹv[SELzSǿ͒#PK<6\ΌTpFͩs@;խ1  Kkoטۂ6Eq=B`a5W]}ugA.|` $ ٝ~auS\-`J Bܨe})syf?_zT{kDleHŪ]7vfpm:u H ?xep5qY1YB4ojuV{3̒A鍣Z^f۸B~#Q!Qe~\I0Go`A߳G^; ~(P WSY׺kz=RI6# +qG'Rjrq|lJve0QHE$Fq D{НIB+Yl+g;lL_< :\솻QsgI n,`o&uEuL)5!scٛ$v[ Ym~(>D:U7\yuDA?܅^R3 FZ@",y!&XVhOqJ(B&9Ca tdم߇]l=8w{i@K$1)++ "ϩdÁӐ$]9ǰÊ\2z.Ȩ$AWWz,~8uW9PQj!Jԉ&w: B0>UmYQ+Na$=tUQĐ/?*iʎnmŒQ'DK\/`c`U`ֈJL#7!xi`*MpflR sc&؂xd!C9Ew]D`0>wJ: +J &L|jUQ=$6W8iVB/AW7J#GlǗNk2} xBL?Yjk`6 ~qL:+<#YS{v.CsT7G&+iK:b`wY[M8F dY.L\mG XsSW2GX|n=hl5FՖ5)lJd14;pS{v-2?uWuCa)j34[P}.XxY0gUN3uX+gu ýs~hgI_T mPv 7[baukQa qγ|2_HU OڻG퓤TQ\Β Yƶ ?VV1oo%Fw{D{Vs}G:#QQL*9IAc좀zɱ~4uODǑә? kw3rGƵoSq{I%64RLh-#76=C 2^Ă7tb E?ёy,,R .a?BpcJ@=i~)¶sv t.W{LSG,xچ1H jGش?6#<F8rb~pI?9,2Z%fZ#o~clc m(fψnL $6%D 27wkNkdmy*]_9(m\}zK'@z!vQ-ׄBmV\EsD-n}By@K2+Aʶ0A jp*ɚmZZK8K,U'8+тXߐ4!C}&/V"[uoQQq_\K9ˎ/'vsٱjV$,yO!,1HPBk:wП@CУVXbc 7GߴB]Xʼnoɥ w1D5bAH["^4d1 jNXɊP쫩.ٹ r=?lsJʽ{7Pu=D]DrZ·rж;6+Z;F,С{@Q툹/ci < W'Nn8 bHHiSD?[>RIqpņyT8 ܾD]"EcGolDX@OcGKwH)8m?b=aYn3Y0h(Ċ2CCfi)*3~;ϱDOô+[˹W3z[6* O'0-i\|!T kSheҤBVcܳ:ך~]8!4yݯMKsG%E0nQ_`1?\S@GH6Ej&Ns@ZVdh^~m,shbt>.[5+ 12K=P5'I`KMש+[>Wf L*_D[4#UO xʝ0}R-W qʯǟz!LQ4י ~:/d׎~ZRG vZ@W*xNI=ԕCk(73 F/xƆI 7&,]~Ir_I}s%?m 1*~9h7! yT~ RezO q{'5)fӁ m:.sih%lHO>]{,F,T?FZD:SA*F ڀB2fk3_Pqt: ŐhDr^S3[ "5,(hmߢkK6#cFTrK= җ.qOX=3%2/kg-5[fNqv2Uǂez{qJh:]D {_qI34m[QfG(|!oUvsY!.5*oo6\~c[9Qo/}!7u*]*VHja݊*/#%(1cE1&~6햻w/"@gڒzLE$~s8GHPsti/ACEZCV5}V0'4RAlrT%:i YS{.lw_'4H>;{<=Y|F wLҲR:٠7W'^~LIaE$fLO`CW,NccM6iF|mο)-.bfpKm(NUT{ÔG7OMWn {CWgZ VٴL-% /__Q˗fKHe_=VMHw+bDd|o32rO6F6k ' 윈PP<[ɻ < *0*J}L앑ڼ: o2S!&t·MVI,nn |,J@]% :҆*w;=g] tΉv:}Y& !ue&f,X-ӗS @1=N $O{δevNlnOp#4hݡq/֬}m }댽ϜKCWwb |o!<=]A* 7H,ɷVC+gi܎/x ^YJͩ\FgMO%Z4( "e6ں{f@ڤ8f'am.0_%FPׁW3oB1 /ё3& Mзԋs+n58@fSE4草"v펭O<6&nsmutNlfXt87KKq8D$_0+PFn HJ%a%J(0Î7P^Iۆ+GԴu~52-0t`)_7d ܷHIw ;߆G\::ht-h-i6ʵؐYoϿ;9x~6q%ILڍ8?,,Q]}'[LĐ>\7֣c"g}exqߊԐ'*(愂7Wy;#gNբxH!x/A__8T+% 'brqer:rr$ +݆|LOX rmv=_9.[Ez  MG|kb׻Pb))$ԧ,yngǑKc$RL3QC~쭐nL咦4?J[ ߧXDi4/c&,-T rJ ӗrI&x%s)!z{VTUڥ9 uId+e#V/LZSi!|(Nc#˔IQRh]+%Dס D/=t?±XЌ+/$ϝhٚIC+h >J!t/eA?Lin [~?_u.Jzq +8'YG9nsu+:" U}"͏l&S0$ 4h2crh Ls?)ZPV2 QGϹHGk.'..U|(D 롾Џb{|a~d!j?=bMPgi2*r]q.oc X)c=f$FlINDž{FjSQLnTj{YwGZs*+tcxߊ]Aϧޔt'_Nf[UA&+ 3VdX#3H2?cM47C[_ aFW+Zp aGfsMRWSӬDCӎ`;iD1!Dɧ5mBr| ^е)&1o+~]x %4`}'CwwO)wR\ߖ ZPl1cJ[5koOM[L*_~lC~iߙ׿_ޛyDbOcU&k N8u E㕎<,Y_}iC_^-2@AA|]Jé?a5ifeO]ܱ)CvwS.Dmx'p'\F8<~6xY[S56S..W B Uۯsʭ»֚qGo uf1fLHƗӮ!dCƃ 3REE4 ߇ѥlgaC~Kŕ>H%P 4U5A>?k uS{AŰ}1m[4*2l bdvy] Akh5C*@BtMJ޺ &8R4xz7t[">Yp`USDW&oٮTs/v)fu I;[uD*NoX棔PvO<f,r rEA:Uo@6ؕr}E¥l)n9#[_% R9)Ĕ" Xӎ[vqF6HbwУxI[ oDa.)6ɞwaA#g_X w"ʋjd ƼBV{eZ uqt*LE>Ka%lGGȃSi]LaE$k]n&^T >*o%w֎O!p9.Xe*{3X>u+/Hi |OW~LR L .A1{8u?Q޹">Rt,ӏK T]GD"Lhy!̃MJ4-Yse1/Q- B@4|%Gv7N`d #%:U=pDԗ] AߵnҫQLn,Rp\rFА44L& 0qO?=J?lCw/7ou|s'*nNh TwDHjoaM f)퓊?J9jLDF5CbP, Vs0\? =ِoHḐ1p)T oq`Eex {9a =sppLc΂ Vɲ%Ϟ}-ʋ YkHV(_X  PG؄8$H$dorBV (WOp#] ,~N)YR4>ǰu8LVhfjCOoG.횩!1$wt,TVuc!n邜LX_c& tK. cg'!}"SwF G>0BUCS.[<2#n~40.eP&jr&Ym**k!)\l4o PPٕw?Wo SՓnfWFM'~F@^jXc8%H|X7z;"%~j%8KPA#} -HF~J$G eY: G:>}-Ȏ"?N|`$lWݻPzNT$E|RkWcSt 6.ǔit77r$5bE.:ie(͓MF"*p2lDWKy%,=}W=9R`O QE`R1h@qky,Baњߔc=7-.̪+ԾhM0{GrSKeL!5>J-INQ`b5W&hqe(R.wS}b 2Йy3aMU9GR6&PFL 11VVL{[~HK`J#W64j1~?7PxeJa`.hua/pn|2סtNDA. ۶y _typ$);<̲w m5P]./)/gC FJOZpw*qBL}b B&E2d&~xs@`pӫrlUuw> ϦS1c"{V(Nh<ϕP!} G f"WA% q:x͝΋͕"7, $Fw2Bo Ԩϐ~Q5C՟Md '7\ko(Uuڤ}GLJa/CR2&u!f#NKuR(c#՗ '-~k,RԌ|*$^w%dz"6pUW>>\Rʻs5A=  (R~gyN:]JhM5oCAN^zyGJ1Gķ,I )dI/E]#O \-<{ЌZX׳4 *pcU2$?u)m]xav]r>Ši2=l'3xNK7!Ĵv[9W<*d\BS0"? yN ~_mQ̓2 7^|6̨T7d_} GtJu)\j6%] IJQ[|:~~Lmnǡq[.þa|}Q. (\U o".uT\K<,F݌D^7^^tFx"Y> endobj 78 0 obj << /Length1 1630 /Length2 5565 /Length3 0 /Length 6392 /Filter /FlateDecode >> stream xڭVeXܶNA`N S@qf@N$nIw9u1^^w '"nVЂB"2] J4;0o,WdB )Gқ>>Dʉq-9Α:*_t?b dRU_p(4k ~x#W'*3Q/ǭ9>=.HIEH ֪4q=0o^ ` jw/{/i])'ڔw;3i׹?Jmy"C1ţ9i*mo'Xc*:7?C&VrjLf q F$HYN<[7 szd ঐxE7G9v M|.}@TZuh{{3Q$ kړE=![la/RH)ܱwta{TPfx_Tja1(Nq- txhƼeΤ**Ѽ_MQel g; x,}jA.q!~>.>U_1l`STJ7?RwIUz\#$~d%+Zx0<*fNqsDyaYm]AdEI}F?'܍Ͼ4|Z|MV>mH!=#)auh}Gm,vT9 *z\{cVMIyO8OQ['IT'`Q24`v"Lٴ;{h#^Dco%@nicۙAwvkhC̋"zcQUX_r tM12Ns7+C5S2O=@ 8*%X\ F3eٴ [qIP[*d'V_My܎0AY}ohK38n+®ϟNz˹Ք4K0"l?&A fgCqB2\h<7{OfKl-]ޅ;hHy[B}VyOKrG 5Yrőa/_24w4Xs尸x"ʬWK3z+Y9[tBxә#Ɍ>_:g 3H83R">yP*whCRfq.m䛿wOrrE8\I+Rg AkP*s}MSuܗ ~DfShȇ9`sS⽧X+ ܔ:&<gI"#UM9vT4_aD24N$s `}MF޴67k 0Z>8D庴T~-P!ˇ܄3\P2R;zdH0?Ifi}uQ]e\/NP @[NXNQg2֮*ؚDlU:nW_:4cK7Uܷʕw_>Xnjmed\USLc#7gPPY}&P:K9u @q6Эw;--m/K C`ISGUm v)xc9| vJ/H\ iZU5N#1>5> ,׋ '65ͧhB nf~P7&iU W)RpVm_g"*(K8۫6ѢmQk/Ip䇘ǜE3pN5JHkG3ܴ)eN=j{k,[ŁgST96NҜA/,15=|QuDcQƃ;撚D{ :Xc @=s~G?k-uF^u./߅8Wu8|D~-J*r5Ғ8{Md^oOJ\.[mre0'ԘO(-}Qzk}]d}sv16'SBZC+68;ahoVABrIC+'{WtcMj6NpU(%$]&/ŔW. s'S# N[Y]'>7s$YKͩƓ ‘ggvSq)mqlU3.mŧ>]Y',r ݅5ϥcEڑo֐^s%I0vDSE3{USבs9wk=IH )yZmmk R:.[@fXb y=D:QQ(#-&4ܲ e8@7}±$tc#ݕu.=!7s肄Z]u{.,3'lEE{[ ?EbrСp>eAnZѻ zы=\s 9i{_ s"KMUndxƉ74cUxPIr/aQkDfbjN32&) '.xջԍ@#W&ϝr'L$ky >]/%#rsUBo3墝D-j<EGlsV"7c2X*QjrAӳ$B򰇩d:YWX) W.<̹fO;o"EZvJ"}?*\˷{/ܺ6O:vJHlz1uN(SnVgEx?6&DsCb\ݒhʀgj9~?:[؍1\g ^oNcP%N#. fZ*w3c24"(!Η:WVfccjG6T[Y;Kf|D fps7DžhI 3Qd- L@8UA`i14 4Oݪx'@M+UQ&bx f;wu2n>8'<=TD>n\9^8a41q {ր6hJ*1*.蚢)ESٵw nCw;D?+ GX6 OYZ%? G^TE;&V_8E(o׾u;B~R#K^1:wtM^^(IbE 8wG+k_Ř}vxỵ.o^u›|i7N>Iu7qoR0 lfxn==dƖYED,s>de8Non$z LDgyA #}j6HWuWLrOV6xo4tikȣcGG4Y5Krם\;idtǀ>a#%ɊGPw?l(uG"ұH|+i3TR Wxz$fX;QHې:h"5z8́X YLw`YuF5 zO`@6ǥ"hpK3j;*{UtfMbv7uN&Jm6\|V؝r]Ⱦ5b/y6N!M'?o.S[#uVD,y,3,,l7Of~iQW^"K񈖄Hm[:R?\mT:"7esӼ_-$|Ѱ$9Ϧ#y:PB:7MiÆщBd!rX|`$eztaTbX[^'l?]gY뛜9~Ao{4m?eK,{'k~!$ ҈:ďL?^gjgR&Nox"?x/Kz#Wvf67&p^w> gRSwA<zFVòRWԸ/'\7_|5U;\Q~LZ]H.P0]XQ/t3u L<\4W/7~pES |c.YṼKR{ތ+___mBЈD =)p,@#h(cMFl$iz'%mm&Fv54߷5U[yݬOZnY endstream endobj 79 0 obj << /Type /FontDescriptor /FontName /VRFQRJ+NimbusMonL-ReguObli /Flags 4 /FontBBox [-61 -237 774 811] /Ascent 625 /CapHeight 557 /Descent -147 /ItalicAngle -12 /StemV 43 /XHeight 426 /CharSet (/a/e/f/hyphen/m/o/p/r/s/t) /FontFile 78 0 R >> endobj 80 0 obj << /Length1 1626 /Length2 12112 /Length3 0 /Length 12950 /Filter /FlateDecode >> stream xڭxUT.wKp Npwwihq-< Nٳg3f\Z_S[MC,j1IA]YJ`;3Wg5OYd鹑ihĝ@.`  2H;2 @rki3|/_.3Z"ڷ7-d\@ - +$ VHANW3[0A 7՚3[.Qg)< LLvX:ڻm]𦷀 af{Kqvq:\oUU$73bi߶4oVS3W-3`kV-o`{!`8,MmAoir5 /ݛ:8z ?1]A,o5.o-Ȭ튬?mn D005zAȬJ,{$/PB ;GO-jkdj1#cjx3_ 1z[b[ɺD6(R` h0}zM{s-#0M @JHJȨ~oߎ*oK?ZQW117%N/Y 뛍/H!ۦ/3[yD^]2]j OJ ;4hWCz>mU fii\!'RaeqxAmCIx,1y wp;J6,kXI嬳<^^wfsj]ot8t4L:b?CQ.^HS)1VS K.<,-$gD? ϯ9I:4 |#9uj=\YA{$6d͕up%?26|+W`MӀ4@D PSbVu FJrՓ`[{>J\SD%|=>(T ]U,`i2IOov$BW/2y Iw3:g"F( .>]FzG#7:j1a&b2Ƞ65n$0(jݛSWX[1΂\Ňە#‹lK<_d5!Gu!ӊ@)h1D~cͫX"nR6g%pA~Aўc.OD|kPd!:B#Ymڊ@qL 3Nst,Ŋx1=r=+<\ #ѫlQX w>6SC3-NZ1xk59WI_>GtNt/^<5+)eLcayt9Y;$z6}$lAE;GcnXQ?dRЊcGk>,VcMD|^k8'L[+EQ@/V9'Uh,oi[;Yjnʛ }еa^Skj cf&啂d9MhӱU9f\ꊶsav!lhj;aAt;ȻbߪDWYGG )z >Խ߬Vr5l _s̒տxɱd@;&Vm#!zMi:~iH :2)qg)Pk( W_V7?n"3ad8\g􄯒v:QP-BlcӪw0$e0Ϭ=[#zš`3WWq@ k2{( jY5i&TJ^矗O .ED^~egp@WRk.ur s4Wh,XyA{-ڤE>]<%Оs'0OM;J\@ܪYsj!<]2Qp|c +#ƔEoO 'P*&P6YW f:2w'"@LU *MFrVyė۩_I?.HzY<k.KYC}ʠrَr0E/(}N> N;肕1FraߑqZ۲?7X6cAw(M~h]X$Vj<ӄT]-&aťJ?l[[cf|s߾h/{aR9iB3ǷP>Tʚ^i8?oP{260-KwA`+0 <̦Vc 2_Lc(Ѷ`(coz_b[dNQ\3{ص.D];|Q8~tꓫ'][Dn'PSvVVK]nDEW.uV#ϗͧҵڴצ6|Ϛ6Ddv:|o>54\] ˈ!Γ'*!^7ex/㊈aCm${JnMʄ o?LpȾƎv e9=WZ|FoҬ] 㖫UcS:N$Z ~9j$*14⎌v[_}Y mAĀis&ctDO!V ·bJ@a1N32DcUZ|mz.$#ZS˂U[<2>k'ȊvӇ8 w7Nm3%&$C o%}o<jsc牣K;ڎ(:`jV6OKfr}Iv s ]d+N+$ߵ5zPԺ?><]KP.}L96`G:bM%&ë'd`~+a*}іȋE8lyT8hF,8g>xψ.90U:TۅҲ%+zù,.-8Z+FpdaȢMF[S)2?%BtEBW!哝aኹQ]K楐81`ubM_}@GW?:YNdv H3?Ɩm_hpye_E9$u+ F\w`ڝ.Ial3 mL 24924RY3\g1@VBƺ;̯5YD驄ĢMbVKmĴؘؑy* g*Dڰ+$ LNreBM@%Aٲ) Hlae #Jo5'A{']Y380 Udͧ jۗ@9X(fX3i@IƓe*O<["uԻeEt~C}r1 PȞIyQ?p|" ]fм釧+3PXO?,TyIд EkڟY +/ED%9yVI՘g,iG/g/:Vt} fNeg(aQ*1dӕ]"f@?Ma?*`%lip!#YJJ0bшG2f+VM*R5 PD[>*&E}kTϟ Hys+VLZ!IV|=Kg*3 8ahY˼ǰaT*r+|vlGJ;C^xA֪9rt2WNLGx-fJ$0jmoB=D| (~bD_u v\h^T*AΛ1rKá̄[&ɨ㾴 }׸L>nb?)ME|`v͡ q_ۺ{J{ zi5F?JD}wJeʷ+u)SAHc.!JI=~*H*ܔ8:D6F98qQީ,SV~t ql*J~QWWܴv^{ƥ)>Ҭ`:i'}3 Y%JwDSE"uuNnQAv@9bI\?YYްb"\aYWP| I$:HNא  >s KzYY-lc|?Y'#hb^ /g #ؔnX74l$V*x+f_~s]rSfa1 C]Cyy08/x%? z/*;񚳟B &uN囈*;2s-mit fA(SSSN EMhzH}v GXёii$p3%8=Je^oe6}A߲plދdK%d햓'EG-|Tˎ]g\Il Nheg]㞵"a'.t7@8rrҤ!^X^[ "x8up-GLH[k1uȤ vPPUˣ(eFpOj@7dђʷhsCH*Փ8gJՠW2#NJZb/h'/jP踲И/T"7(;I]`MGM>`RO`㉸|Y[cJó>x Zr%C!e. WT=y=,\Zu[/"roj Us"\b>b*/+Ys d@/[߬ډwnHҽɜْtN4Sɐ;O/w<hԒǀS txPiKz["%" 1s bᏘ #6hBW1  g}v sU'JNLR[ w9^a$r"-6 *+\Î&7Mbr|HUrs-+tCgjM({aL08zٍBf0/<\jո$ࠑ80!vn~~]V1W#cx7EaSVz꼠F.u#J~ǯQE^dtkW.Dy_:3G?%[p}(946@moSF|zG[2RnťA[m8?ϥnh +K#YsX&1?_&hꗒي =]^.E%ϖ́>՜ 5(Jp)h:,]Z~ hvA7W9<4N;ic?{QsN9g9IRS樝b?]01%(f3\9dNW`D?C"L`ȽMP=D<&@QWP;Smf5v_7i֕"Q#!.g(63-0-ƝP.|W&Cƶ[} Rz5o%.R}Q=¤tƕG8sj~jq Od)Ǚπe r?-86 knfxxY)+!$$ {_bdx2.1>-_3կXhQ~?`e>]@piϙێdiy%}8 jA mH*1t |9$X9M- w9ٳ!.rgE+ۧ_{{XiːjDӄ@: }ESyNmV68̬tcjarHgX?ȠAwf'?6 ̈́D`/?WLtaB]:M{fPbޕ$S`p4%kpdbsgMV?~>mVy UJ*GTIlIډgk~0[4ɣi2Cx2GmOҁNLlk̚,Gܕ}9w HBXfj'IΙYGs_ae;>D{5{`Anb)!OLRZ#A473:z}SO?LcDu؛h̩>~lb#㧇ꩉvd=J UE҉#pnjڨ#8>kw^ϼv)9kB; iMfEm/?c Hd;CnTƜa@reCS+&gR* 7;jـbyЄe 3E0j5;+.ʱ&!0h o+3UdR:~}iEүQq#S`% 5r}.qhGr N|(s f7t;>oO"40uv*{նct%~)ʋi4-l}yN 5en'2T]I0-&W*+'~}/Ϗ^[%z>'S77}7;?p&bVQUHP$f%gEY(hbMZ8S 7q!:T{T!n]aJ xߌH&hu''nLdB$ߏs2=fT޿r.*NՎ+mK3S? m&܎tG^h[B9҅V$;'0}:ɷZ}Ҟ?s\|k1D#J.Ǣ\ S9FSˇlmQߎ4{w%جiCK,N]cd?e_~;)EnSp"\ s}/3)8Y 3sb;90&IQ0yִղR+WM8ͶTDzFG|P~ :f-U".].?l=iRh98+wF^`mk1*?r|񻀱'NCr_tھX~G-8kJܹQWϵC+KL`c,MژCE;_^С` ?]Ъ!lrk/킾U3Xʌ-?>PurDe)ivqdU3[1u/S .% *.fd*?vCĻFb!cؽ̟iI4)UOg+a88`n3Br*qX8%O&܎}epz?,:`mJ f$E/@R9pnWo@(=O s;XvgƎ6=9.gCK%jsgW@ ^cPKiXj}/e|BYQg+O[0- z<ƧooVN0g4'U^$,D`npd: QX^uF / ͰJ(Dd/pzn¯dA<1_׵aFh8 )M}d B'Q?#Zk W+Д@O\A6c!\}rH6o}L}dQ>p@0zMN|ݠ'<0ɚ7B&{d)XNgve"uʚw@⒏V޼K-d'Pg-BY>EdmeLf":'9t5V ݀@}V7Y`Np,CMqޏ/Kj(0Ky[雫W endstream endobj 81 0 obj << /Type /FontDescriptor /FontName /EXEDHQ+NimbusRomNo9L-Medi /Flags 4 /FontBBox [-168 -341 1000 960] /Ascent 690 /CapHeight 690 /Descent -209 /ItalicAngle 0 /StemV 140 /XHeight 461 /CharSet (/A/C/D/E/I/O/P/R/T/U/a/b/c/colon/d/e/f/fi/five/four/g/h/hyphen/i/j/k/l/m/n/o/one/p/period/q/quoteright/r/s/six/t/three/two/u/v/w/x) /FontFile 80 0 R >> endobj 82 0 obj << /Length1 1630 /Length2 16226 /Length3 0 /Length 17077 /Filter /FlateDecode >> stream xڬeT% NpPwwwww'Kݝo>}z1es\*JR%US;c3 #3/@dcbg`#Ǡ4wsSR9AvF@^& 4Xxxx)bv s g& jqo@hmgou  -35 - RPHmF%ck @du&vZsb%08M@Àn&@L{ 70w4u;g;U_`JvNN& {g߬JNf_OS;Z/_ ts'1` r6r/#_e8lz#/?>KF6cgaons-(o!7r+G#}ZZ#[;sh1F6 kSؚ߱J$ hr6Yٿ궦@Gk-/+,@&Vov$.')#!I9rTjkۙOQQ;7' '7[+&忀XS7rvt?%#akbg:ٚݶl]݀&L,2ӝpr'u{Y KՊ jz|¶x* jy?ݗNdhG{/ i ~Sur0"jFy^.dߞTV/y!dsrHGZI`c *P5t]N,%+UCzg ~HހӾ*.X"_Px:^o_iHiOY3[Qu4W˧b 6IxgDTغ ! @Cb\?f0W:67L']T柬 1YlhNh;& c]:>v5H$-4"ZL^1cgI$VDEF/g1UBAAJ8'b1;) gIVߠ84OW1sx}@tOȬ%栖;[oK,Bۋ೭>>I _W/m #_,bP72F6*G(< AF5+9] ܄32Ys͋ R5N3;Zdeط(Ol*1skL((o-au,PpYE`GHtNY=$Ͻړhm[-Y#_ދVcc;RSTxB+CP Ge>_7^ŀ]?i9kNZ":RkiShHqA글H=![ /yeSuh out:ѽ֐ۃkn;v|.;uF+C0P 3 \飴GSH!; ٯyfW΋x}Xe «?Yﴋj?Q/o. cCQ1 ix7~8) J=u>&c/'Slg@V^A]ꅟr3dWrH,q.A`J8æ׉}tZZSqGfǠ39$³BY8Bj)h[9)/ޤ5CU$r(V; SvbKO6bBIGP{B6AMc< 팗  Fk]Ҕ qtg樾UgA׬aͰ%K0Mr\&4YLpmʎ gEkِby1'Fe/KNb3f]}7sZ+?joPQG|lM}n-ӓ$_%rA5n^{Q8 OH<t"Q ױx.u18$a𲌇c#N]X؟S6f'c'cf{;3Ԝ'($ t=V+]wڝ+}3*6;ſ ݛSg$/EI%#<\1$FZ\ ŕط@nPsLye D*tRQ3J%5YOl]#i(( %*K0tC!n 0= Z@ .a&e}W?T۫QQ mY|+w/n&Jن遵 HgfeظSϏlS ;)9"iTY-fu 6p-Ir$)"2"ļԀ+#Gɼ 6@ݶw@yf7 9}-ҏy\d:Mh:'y?FlA(= 8Y?K[{*6*>y:ON"E(1WLm'[:}aƴa7+ϰ7b;pv'1 7-a0Dfa{>Ad쟨:(LL}%uH1Gc{fgi=TsJ:,3JC9Z3};*Wq*atUqg'KohT5 *YˑWE'99~] p?ea"Dk,4`ڦ#yJuv2^j0o/Hi׬lg:*<"@nAZw!2f,57\$L~%z]"T}ooVv+ /zg쒭 1v9̟,=c"""J})ǩ_/FX-eCu_?H^0>=r XF뺽(i^F;ŌZ䪪2EKHKAvH˹Fѣט ^JCĄ_R#=xh~&]!Iֳٜ ɢm et(=}7MNJCFV RA%m)/#)|4sxJΥ4UxXeyN:Ew91XfԪ8Cۘs{j_ӆ CO+7/;o2ࣘ8tu[EGͱL4۹:K=/h;q)T؆?<rt@-E6>呢}{TKi%yLl(Q.B&Q''3.4BFR4u<.s3dS \RYQzࡉes>^!Vq\qj;im~AU45tX7ͷ!VI49*ҡŏ;e[\4U E'ZYoǥ뵏~&gq:%ҞU繶쪕Ekf`f7[tBAmNGvoQYf ]T 6h&VNyWa CFBAa*Kрg`hC7`6biu{wq+!eSƓ P=% F^@ QVyWxQʸ*8)Qpu 5dVJZ9k $дpOf^V3G(In^@MJcFӚ8V" | ?#-\\P F85$=ߦyRZce yChEX䐁moz)Lf5}I=*> ,WHtiM.e(N]qHRcq.NSQۚ lz;Hgx!_pUM^1ٮ0* eI ?&[.VP1=s :{6uez"nxX֘/g8}k2o.+&eǹ+pboI?o &3}&]Y~lXMhǫ79tJn#%~Ԅc2r5O 1pUf戓ijGG,YL|35n(X'B+?<euP,@q+ȩ ~!L4<m]cLT(V؛>B)"n5yrJZCZjYgpeBP\aAe~f{'#QgDWRǶ\n1m2XswG*Sx#&$VWh\[-tT]36eXCGG +Ziiim *@?.4&̃Qd@7{ai%bd5ݪԪS$Pg;-d ZbXyӢL|OIT\@Z{|JO_|[J|}hl'K!ۘ ~C&+GL^+RWa@2mGxb:w)FGk\c+$hV9}J?قm33X0%bau rMd?9juq9 .^]֚P(6hL"EwU%VJ!-hȣ; ;§crD6~)v os-n|f CtŌb]}VmrRY{>@8m@C)sd\ :t69^*=ΊDLEx:TDۧC+76 uyA}zGK<_ۍ\f |YǚѥlD ،LCHpCd 3ekW̊_'Nh1! ׷Dk*#HkC_mǧ~Fې4Jߵpoy ԭ/Q9 QJn7,+̓o ;,;5̧ @_+ т[D3C]^6B]Tz66W.DUga`C=L#g;Ո<z씂TȏȻT>fDt'me۠mm2(z TLG|2wt X8[Lnȟ#t\ -C_7|KUKD>"Pw)/.}W.Iq:JU!G!#Ct3;͔CB#nz$w5Pzqͦİt8eyZNM-$-\4*b/>lؚ)ʐ)<889D6"B-80P$ :4٫A)>a|97꾠m7G-⎇6g.7 QM VHg}8muPu녥{I9.ޖ9G\G)_H,xIne !Ħ^Toг:K*$~էV畺2\H#jSWD#7)ZmLS~m@c9]Y4[M;6f=7a++7]#x3 QD v6YV76$u_Vʲ"B'(';*$6(~ᆑFTr$5NVn}LZC{]ylT Mm4.,"Y5, l},W)-y:Zwڰ ~qLBGN,c :ƧmWddz"V|b(gY)e:v> _F2p}> 4仕.]k Y-=ŸGbm}a},pfx cK Ar,KQX$ƉA5@rzGڔ M*~=Dށ䫰診<<`z2eh?{mpF$[lAk 2E` GD O% ^$ &1 N))jk h\'(/ELC7im`=ʲZv'6`!ydvZGdCܕv^-Gx ׎Cv` ug}bŷ@DKv)|<~dDpOf!&TMN 2xAy6`d `P9׃ <4vd=<;T#} R]ѩZH <13-sA"YONnQK-Q#5$lHO8yF6 XQNsV8݀CEd$ғ;4;Bq+!%6XXS7M$}Cs|sQmOؿKTlrT.mJQ)AS'eЎ1Ȉdh_,Hߖҕ(d@S/[\s6w+2GS LΦb 'yr>;$+k4 4UtBsd,O;&Sy@BzEb`zs䒜Ly9>%suYfӨ|>M΃&3>jѣ%BnqFFicJ Q&_g?0eG4-k|z uҌ\Ԛ-OS1:ދ'< sM)T̖IREJw7{Fkr@mF%Mޫdy/FЂɂo: BBTxZQ= G\}fy9 ,zyw6}dqL˫9}cD%>D~9*_+!xTO>LWD*[/}ң71| !t8gEMv ۬0ߓm-)ƥi{3߽<唾Q}!Ԋ,=.So^:GV Xl4t=2?G"P*d=ugsH^[B1Lۍ<$1YP!!Lz.kxO( xd"eoru<oҸ Zn}g;o_v;4a޿5,8||Oe;2_vI>0w1sKBwAe_~LcּgxyAtT@O:6(tV2{ Pk`OTitWǩW6ϻ3d<;tl>VEK-vI Iqi 淏VP*X,%v &'z"&;s ZD{m;-Ԇ:AsfDOܑ@7b|JnUQ,> LQAl3͡/$w_]m$c kI/کBOEI`+"fwwp-$k]mjy#JHϐE\^ )Ic*vߡ2Mm]Ev*U3Ƽd1.CJ,:FU-@} B,Fڬ*ǝ&p%l벎uH9Х6EΠ&~z !"DK[rqdPb<^U`kBB##n`GZی"֜#gM޹lGR,n_1|+30[ FqINx5fY?j a_wZ뢸 ~ݰfaJ)>/b]3pC&{& +I C>.0Fj 0$m534}~Ap(qljΏ9bys62fhkn.OEĨ<-ߜȪӀ)283ഷϠ{0˖*žx_XΤJm"DFZQB[ߓyqҝŞ> u.E:;s_@4}.uDzM? [" GM`8=1ۖs*VjEyw)@[dh x ? KT1SA5fsu]>w%%5v2%v)%uRf^I㭽GGQ[e*v꿧mQW$/h_@Zǝgza)Ā E7\.U4 0_P}m糴 N&~7Bɵ_q}@sr0~Kya햊9qWҢjDuNyͲ!'",J|+2x^ϕ4J,PF:O0v Ne)'iњR!&ma$}ӌYɑd,Ow)8pWZ(ca%Ch+`L`Xo/OȞ ߢt.Y-fW+;9_ 6(MhtMh`L98U]'SpHo`o5|=jD23- |ZF|Lђ1MԟpSe453V"`{,8CPM 5z"\E(m<-bN/#wڍ;=Pl8'.Y St6.)r~̏>BCÒU[m$^H:Bn\4 O.kY =R3Q|2lbD`{E|Ѫ[u=?zY\B1d8cDɫ;,!4L>]ia(x]9aWj)V^-iݠn:Zp{ԉꞵՄ_*C[6<'s o[9 wf>b]nޭ8+τ2S~ 7AfSAɜMv&Xt/hsv0IyRvm@-DVzDi(O0氥&Q:آ =ڥ1![c|f89#]~&wGl f\.+QXxLC B$5Aݧ8b7!(8G`lYЄߦ]=.#{!H^Vx[_XZ3.7b}yby.BpILE[׽VW1sC BwNu8ϰ锡Luv6EGu ϩ9rp%lzX-$Ҷ86&/b lu\3z?)]-S&MD_,yS nqJʘ =kӍ (+ M<*| N}fۗ^PGi 5.~:d$B,350o`V(#%CJ )fb-Mߵʙ+g!D29j}U}L“zg(mB& L1|AZw3 孡9$ 8Q)ƠCKB-^/fȏX:31 >xieUCRH(+-)%J,8W3Ӟty.ے'_Kl5?(ai,yT`TjfurHȵd#úFyU!WY=4\X04?hF.34e8jυ:RSϴ[Ikib9թ8nV4fhNoXa4c+:q;` 1E=!'xג}Uے<8wnS}nz䪱D!:Dw?l֕,+<#)WYA7է q~;\Y\Lv94ևaҜTb0& [8Lsĺ@\Em !1Z%]c熿UjK j(vH O /%ViW`SZjVLU;\/zv2ֵNhg{q%y=>{!TA:5pSObiPAlOZgrؓ& `,.e;׫?_C941hRj-j/G7.B鏛uԒDZ);:07&K6PE~Ib)IDY`(TyuPm|N;~[uB;4]~j5FQj{3ԏ=M4x ck-D -}^S JWs>CiUI7T:?`FX^z+{*jWP *hBB N3X]e?vs[j`3Yכ;xwƗŸu(g'"'&j JZy/Sf}o)ntލ yao-<Aǰ$ds<ŞPLIaiȏ *E3{_dM#b|cğoѰ(@|p A&^Q^ELO:Mip˧cȌ22ݖ%(TzhbR8̰#Ņ -%_5}g,3Ze5gM?L<Xct] .Ykz 1̪ItJ.t&[?p9eY [Ț# KMtLn\F;Ѿ+an#6X!<L(I\}Er՘g&|#VR E֪ޔtiD\\cDb7RDv}O|ԣJ.%AdFIl͇f7/^@Bj֍TQ =0σ4||-A4,<"=}DnʹO6oz&G,6!i$:~.9Ņzyȿ(G2fH+7Zٝ3uj yM#BUൠOP WߙOUF;\ 5yj_DSػtg`b=%0) ,ّbe3AQ ~v|gbwUEB[bjA%Ηz֌$xRL.o-{p/iX_s_3@!*ʋEBT쒥h9x#5U@) /͍b+m ;$~d$Yh endstream endobj 83 0 obj << /Type /FontDescriptor /FontName /DLFJEF+NimbusRomNo9L-Regu /Flags 4 /FontBBox [-168 -281 1000 924] /Ascent 678 /CapHeight 651 /Descent -216 /ItalicAngle 0 /StemV 85 /XHeight 450 /CharSet (/A/B/C/D/E/F/G/H/I/L/M/N/O/P/R/S/T/U/V/W/a/b/c/colon/comma/d/e/eight/endash/f/fi/five/four/g/h/hyphen/i/j/k/l/m/n/nine/o/one/p/parenleft/parenright/period/q/quoteright/r/s/semicolon/seven/six/slash/t/three/two/u/v/w/x/y/z/zero) /FontFile 82 0 R >> endobj 84 0 obj << /Length1 1647 /Length2 9861 /Length3 0 /Length 10711 /Filter /FlateDecode >> stream xڭveX\ݒ5ܡq$CpHK ;AÐ;w}g89֮UU{ծSjY@̀0K 3svR)CՁVrPS[ ƇFO/4 `IS(PH<<n!!!4zde 0i0g ȋ `xpB` Q@K- ',`Q@Ǘ"TlAE9dXB/OiN/\bNS=t3@G;7r4C_ x[BJ{!S8AAPKTUIZBvN ^h^P)ݠb '{[S/dpvheha trzys:ߪ7u׮uZrq4Ķ8nl1$L- `[wS} `ߩ E7qUv_mmM^9x4`ˬ( G2\M@r+w1ՋB||AN 7*jn 4}9Z` - |}qLdnl5WZ:ôkKW@5큀*qq_#r_'_D\+BAn...n?WB#6X# )6wvt|QiR?] h6b&cjF(g`TҠ+7@}qfA_%75lM*~\ẃӖ䭜`Vwԇԣ_Ȩ= AuJ~52H6m,0xkeIg=Һ?TA{YPݶQVǵkc_! "V0>a0L|høwcoyݭBm+ŐO#&^+Z g}"~xv;]{T1?+It[|)k/+'\IiA>%M<~2;qs6mGd& pqq|t|Ml8gu,_moeVZPsz_J-jozܫt鈍[]dV >ZÏ&R".5ltRЫ N#) _#Z'\BW vGozp\3^sL۷PX` 4qض.FPvg|Qi+R}80 !ОQ^DH˛TA)V2qAE?KwFw^eVxhÀ(; w0\vZ,G2mz Yϴjd-4"ABO+r1^/;}Aqo_/f— 1.(ǚv8=k]LWp%ho?k[ Tۡ}wiuy7M0yDhWU|&B)vȎXLRht4QKra'OX혳!Nr=K Ef_"s!RBӂ0ji&%'넦ֺ8|+tȇ)6oMXK'2,E;8m[Jե[ûFS\ָ g06s8Ң"!7o)1,JأȐп@ƞ,2ybεz#c{/22;jF5QedXCR38{&,IK/?JGIJ^6.3 {%wOeH}YΔnv-$鼯$H;A\p\dB:&?܉ꮝH )(nFUPY7Ԯ ׅͅ "LJ٨aԛs붮2kChWM> Oe hбTs(+!o]&~^>-/쿝y R([pQ0cܖz& m<-% 4DBiPsRʗzCm-:(GnKS6hQ UCz}1LHDXx oX|$\5-o}aL bh@b*W!eZWg<{&69: +_Fx5ޓ?T\i]8< q>x ^6NZz0ֱHd+;a6ĶU"نP6apvnj{חHx`y/| o>ă?ﴆ.bqkitgG?@u!dDVXָt,taKT^]*vwةO'A2^(RPp@wiBrVZ*@}Jt-k M2;(ţv,:aL3euC`p]jk+̎Qɷ.SaHz+LsIjk3$G9K7^b[RhvYXsSjtG7o8_)\Vz*f@GI2+BQA1:$nQ(Lod[oA2'"R- 1i{d%X%jEB[RTNSȌ-ʴw*nxg<~i.#-3Tbt*Uc "-X(hbV) cP* k 9hm!}&JbGf3 C2xtd.齒"I:Ψ) L̽6ǖ$mV5CtL9yU'i 2ms_9Xjt]+zak#咻 xnLwM̭qZ|<}'YzO%ms[RkyZ'2H_RX?baiM.F*Se&R2ovܙM?4 춆l&3pu5 >& |ˆ]|;XGn!'Pq`fi^{ _"rwޙ33G!4B1yX9spZӐQϵ$ kpV| Sxl*p_z)IA3s&VFϳ]t]LHh|y:(Qs80m07bա?QH[|ꕹ-J7޸)\_!͖C ʠXX,s7$taO.'F|mb\PL[՛q?}.AGaY{v(}.?>2,^)[N& >cKZGBp 2O[V [ 9q Nx{!MRxTлΏ ^.;NOG,xtnTĆ<%Mt>7ӫ#E 6RmC*ʍ*@HoG;31~+Ap ?1BHQ*2r\$i|ۤ5.#2o0ֱ؜4JmW^ŜE_5U GℌOfXx+ͪe{Q" O)O!ҷEe')!r-U- {좑S@b"c$4}AG֕lهp*uܶñ@{!]n,mD8CG^cBbR_@$JBzx{&:;933PW }t8 *k=TGRհ\_DJU{H7I}Xt`6n#L'BH}ѷSܦUP``bw/U~-+;o(;*(Jm^ECZ.~ay+k9/lbO>ˆ*I4!sSB/Ep&+&Vn$Z5'P*Cd9Ȼ"T!Gݗ lh0?JQJoIm(qr䜔y.'!LفftK-IL\=< ]ki)- WwQ1*zKsg[ E7K,2|'|wKv9\XqJA"~^u^(n﯆=lʷbi׍v;EݎW+4C^߂KUS2}6gg/g kLJ( J[\`ۆԋ!JO~$iD|z%Ƿi/]BϷ'e#bOkA{_A% E1<)?syIJI:[5ōQc%l)e#Yʄyh\s aޠtg yh"C2c|  OIҽ t `m'b;Y#fbsx 3:N{1~-f'@T#lDKycw#\bc1M¿J+Wm 0(FԻ#Md k;zZd Sq2t ~̟!O%<#jqvBh&ܔ-.SAi;zi[GN/k| n&'qaL3s"K* w 7gUr̍۫W=C:?>lVyL$b:HS?c}m5pGLY.%DʋX+\H-HR{4Iϋ~6!@߻4`oca;d:r,{Q8R~{<)\ `\ a4pnfr8:vɍ4Wm^Q?6W%k1|bޣKY?鴎S7ʿ3z؆gUH â/[L>֑Ejs=\q_ҝ=9w)Tk]JȆzpÂrO"9'GgL"zyj櫓LA,o+ nsm5L Lm t L-O;"AsH7ȍ`6;F ɔiE1`Fw@j9G&6B^;$G3˜)GUJQxޘyAN(o`:vmuXn(#ůJ>мLsW ye,pj2H[qބ)M!]M f&!]&jx8 ;g5jJs@2-Z}YչH.IAtbXC0=R.p.^ wOoGDd\xlMA%#0ц`a`1LC߳*zۃ#+%زXTm\CEsqA.BPKSbK|=[o\c\R0h1'H=h܌fVq_37hm3R'lhAt7ya^tZD΋2~> 5XC_P[@ayr;G5.J$f+i43N=&;)Wb-hz)2!Nݦy/"w3>F%̑InՀĞI"8cJl`g}j>4/ȃ2[ ׀Ɛp:jϐ܂!d|ܮ|8APg(]a7ŅM$[ ʖUawZyˊrCFSFo 遇@ M_l#ښrF[0R¯DBė8ę7:c';w]Nz bLrYTku[mBF-yE| WeK}NJh%TV9AWK{k5z9p1٭6CcsE-Vz|H6PzQq{H/0iCԕ϶JIdqQlgZg87B, l۪HΌol?/"P^{ZzcVtfUt%[h`JB-`?`E8W⡯z뼺zx턅oLN~c6Ḥ]ZI>fߋzJځDbޤ}: 7c-߸1u!5c?ba-sTr⚋br+gת XDEuޤhTi&#_ቘϚ4 ɚ+"MNeb . Eild^Ok>׍ gVNtOR vRƸEߺmIN("b:_M uSuv,Uҫ^8[ytvd&hjЫz0yrR$;Y©ދVu&|D%̿4`tȨtE6F<Y{q?:@tksqOWdj,έl jO*C|A*`Ţ;W}S*E ص0w|v^ZH6AgVRz* >]G!71 qyΡHi--~fކQ>*>*s޴ě./Zbqx=5a&HnB懓p|)2X5R:.]cu}c;}E` CS@_k0 cR?$(E0ߴV\|Q,vc>FSpUf;mH _64ɭkep}9Xc>MN"¬ˑrR{+wHe}B5!"[D(U[:6YGɞAw wn!('Sny9aQ"X+$]\|\Mr6%ý"&_ed4$v8wRs(P%v 0VhjzC0ر6~Rt-Ͻ.}.\Us:_ B-YJ\l,$zή)9tCnu>Nn*g*9i,U6b7cńVR/!)T~3 )22Fw቏F7ʤsWr*,:gբVbdzAO_E6Ӳo}`YJ5xSȿL0eRn0V,kz D<~,dȔSÑF468^XlkShW#_MT"3s> :/Y1Ajws"do~ilz*o~/4Lʝ|rVVt3[K|) FhkP)Մt*8W|32pu=ӭ!u=ZgEgɆ'D?q)MhUD}_]Av)+`N%>X,dkyߞ+Drou^Ne>*dQ&ǔڕ+Þ hL&!] i zp߳$u$nKN}/B@ UYmepg{Ǵ_Nj1u%q~Ϡ"P&}c`+erH߿HuQZWhI.Xv-j#Ĵ65\!ء 55Zb -NԔNeI3"`*UmF,m^--ixT][cWe{Ûˬ%l8wU^9' SSQ=#0'vXW!N=a9Zn,ӟ]&:`FK`>>o1 z;h(L}]gI+XP++8&قﯻ~&$|jh Bo endstream endobj 85 0 obj << /Type /FontDescriptor /FontName /FUHDVW+NimbusRomNo9L-ReguItal /Flags 4 /FontBBox [-169 -270 1010 924] /Ascent 668 /CapHeight 668 /Descent -193 /ItalicAngle -15 /StemV 78 /XHeight 441 /CharSet (/A/C/G/I/J/R/S/a/b/c/d/e/f/fi/g/h/hyphen/i/j/l/m/n/o/p/r/s/t/u/x) /FontFile 84 0 R >> endobj 58 0 obj << /Type /Encoding /Differences [2/fi 34/quotedbl/numbersign 39/quoteright/parenleft/parenright 44/comma/hyphen/period/slash/zero/one/two/three/four/five/six/seven/eight/nine/colon/semicolon 61/equal 64/at/A/B/C/D/E/F/G/H/I/J/K/L/M/N/O/P 82/R/S/T/U/V/W/X 97/a/b/c/d/e/f/g/h/i/j/k/l/m/n/o/p/q/r/s/t/u/v/w/x/y/z 150/endash] >> endobj 21 0 obj << /Type /Font /Subtype /Type1 /BaseFont /FGDBER+CMMI12 /FontDescriptor 69 0 R /FirstChar 58 /LastChar 106 /Widths 60 0 R >> endobj 20 0 obj << /Type /Font /Subtype /Type1 /BaseFont /CVVVFW+CMR12 /FontDescriptor 71 0 R /FirstChar 40 /LastChar 93 /Widths 61 0 R >> endobj 7 0 obj << /Type /Font /Subtype /Type1 /BaseFont /NWIAKC+CMSY10 /FontDescriptor 73 0 R /FirstChar 0 /LastChar 24 /Widths 64 0 R >> endobj 8 0 obj << /Type /Font /Subtype /Type1 /BaseFont /SDVCAT+CMSY7 /FontDescriptor 75 0 R /FirstChar 3 /LastChar 3 /Widths 63 0 R >> endobj 4 0 obj << /Type /Font /Subtype /Type1 /BaseFont /CYHTSX+NimbusMonL-Regu /FontDescriptor 77 0 R /FirstChar 34 /LastChar 122 /Widths 67 0 R /Encoding 58 0 R >> endobj 19 0 obj << /Type /Font /Subtype /Type1 /BaseFont /VRFQRJ+NimbusMonL-ReguObli /FontDescriptor 79 0 R /FirstChar 45 /LastChar 116 /Widths 62 0 R /Encoding 58 0 R >> endobj 5 0 obj << /Type /Font /Subtype /Type1 /BaseFont /EXEDHQ+NimbusRomNo9L-Medi /FontDescriptor 81 0 R /FirstChar 2 /LastChar 120 /Widths 66 0 R /Encoding 58 0 R >> endobj 6 0 obj << /Type /Font /Subtype /Type1 /BaseFont /DLFJEF+NimbusRomNo9L-Regu /FontDescriptor 83 0 R /FirstChar 2 /LastChar 150 /Widths 65 0 R /Encoding 58 0 R >> endobj 22 0 obj << /Type /Font /Subtype /Type1 /BaseFont /FUHDVW+NimbusRomNo9L-ReguItal /FontDescriptor 85 0 R /FirstChar 2 /LastChar 120 /Widths 59 0 R /Encoding 58 0 R >> endobj 9 0 obj << /Type /Pages /Count 6 /Parent 86 0 R /Kids [2 0 R 11 0 R 14 0 R 17 0 R 24 0 R 27 0 R] >> endobj 32 0 obj << /Type /Pages /Count 6 /Parent 86 0 R /Kids [30 0 R 34 0 R 37 0 R 40 0 R 43 0 R 46 0 R] >> endobj 51 0 obj << /Type /Pages /Count 3 /Parent 86 0 R /Kids [49 0 R 53 0 R 56 0 R] >> endobj 86 0 obj << /Type /Pages /Count 15 /Kids [9 0 R 32 0 R 51 0 R] >> endobj 87 0 obj << /Type /Catalog /Pages 86 0 R >> endobj 88 0 obj << /Producer (pdfTeX-1.40.3) /Creator (TeX) /CreationDate (D:20100409161158-05'00') /ModDate (D:20100409161158-05'00') /Trapped /False /PTEX.Fullbanner (This is pdfTeX using libpoppler, Version 3.141592-1.40.3-2.2 (Web2C 7.5.6) kpathsea version 3.5.6) >> endobj xref 0 89 0000000000 65535 f 0000000924 00000 n 0000000820 00000 n 0000000015 00000 n 0000107179 00000 n 0000107516 00000 n 0000107684 00000 n 0000106905 00000 n 0000107043 00000 n 0000108025 00000 n 0000002002 00000 n 0000001895 00000 n 0000001036 00000 n 0000004487 00000 n 0000004380 00000 n 0000002093 00000 n 0000006541 00000 n 0000006434 00000 n 0000004589 00000 n 0000107345 00000 n 0000106766 00000 n 0000106625 00000 n 0000107852 00000 n 0000010307 00000 n 0000010200 00000 n 0000006691 00000 n 0000012313 00000 n 0000012206 00000 n 0000010434 00000 n 0000013505 00000 n 0000013397 00000 n 0000012440 00000 n 0000108132 00000 n 0000014623 00000 n 0000014515 00000 n 0000013609 00000 n 0000015805 00000 n 0000015697 00000 n 0000014727 00000 n 0000017647 00000 n 0000017539 00000 n 0000015908 00000 n 0000019241 00000 n 0000019133 00000 n 0000017785 00000 n 0000020610 00000 n 0000020502 00000 n 0000019379 00000 n 0000023081 00000 n 0000022973 00000 n 0000020737 00000 n 0000108241 00000 n 0000024236 00000 n 0000024128 00000 n 0000023208 00000 n 0000025701 00000 n 0000025593 00000 n 0000024362 00000 n 0000106284 00000 n 0000025839 00000 n 0000026305 00000 n 0000026601 00000 n 0000026908 00000 n 0000027214 00000 n 0000027238 00000 n 0000027397 00000 n 0000027967 00000 n 0000028435 00000 n 0000028809 00000 n 0000032182 00000 n 0000032427 00000 n 0000036289 00000 n 0000036605 00000 n 0000038257 00000 n 0000038509 00000 n 0000039873 00000 n 0000040103 00000 n 0000056863 00000 n 0000057307 00000 n 0000063818 00000 n 0000064073 00000 n 0000077143 00000 n 0000077503 00000 n 0000094700 00000 n 0000095155 00000 n 0000105985 00000 n 0000108329 00000 n 0000108402 00000 n 0000108453 00000 n trailer << /Size 89 /Root 87 0 R /Info 88 0 R /ID [<3CF6A9B9E5F2EF62566082D67CE51DF0> <3CF6A9B9E5F2EF62566082D67CE51DF0>] >> startxref 108724 %%EOF CARP-1.0/array.h0000644000076400007640000002623211356774460013332 0ustar volmelnvolmeln /* Author: Rouben Rostamian, June 2000 Revised Jul 2000 Revised Sep 2002 Revised Jun 2003 - macros don't call exit anymore changed macro names to all-caps Revised Aug 2003 - changed reserved names to all caps Revised Dec 2003 - changed array index types to size_t (were int) - added an "out of memory" message, printed to stderr - most FREE* macros now come before MAKE* macros, for possible improved efficiency in preprocessing */ /* CARP: Clustering Algorithms' Referee Package Copyright (C) 2010 Volodymyr Melnykov and Ranjan Maitra This program is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version. This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details. You should have received a copy of the GNU General Public License along with this program. If not, see . Authors' contact information: Volodymyr Melnykov Ranjan Maitra Volodymyr.Melnykov@ndsu.edu maitra@iastate.edu Department of Statistics Department of Statistics North Dakota State University Iowa State University Fargo, ND 58102 Ames, IA 50011 */ /* array.h o This file defines the following macros: MAKE_1ARRAY(a,n) make a 1D array of length "n" MAKE_2ARRAY(a,m,n) make a 2D array of dimensions "m x n" MAKE_3ARRAY(a,l,m,n) make a 3D array of dimensions "l x m x n" MAKE_4ARRAY(a,k,l,m,n) make a 4D array of dimensions "k x l x m x n" FREE_1ARRAY(a) free memory allocated by MAKE_1ARRAY() FREE_2ARRAY(a) free memory allocated by MAKE_2ARRAY() FREE_3ARRAY(a) free memory allocated by MAKE_3ARRAY() FREE_4ARRAY(a) free memory allocated by MAKE_4ARRAY() o Additionally, it defines the following convenience macros as synonyms: MAKE_VECTOR(a,n) same as MAKE_1ARRAY(a,n) FREE_VECTOR(a) same as FREE_1ARRAY(a) MAKE_MATRIX(a,m,n) same as MAKE_2ARRAY(a,m,n) FREE_MATRIX(a) same as FREE_2ARRAY(a) o Additionally, it declares and uses the identifiers: ARRAY_H2RESERVED ARRAY_H3RESERVED ARRAY_H4RESERVED within local blocks within the macros. THESE IDENTIFIERS ARE RESERVED. The user should not use identifiers with this names within the scope of these macros. o vector/matrix/array elements can be of _any_ type. o If malloc() fails during the execution of any of the MAKE_* macros: . An "out of memory" message is printed to stderr. . The macro's first argument is set to NULL. o After a call to any of the FREE_*() macros, the macro's argument is set to NULL. o The FREE_*() macros can be applied to previously FREE_*()-ed arrays with no ill effect. o Note that macro arguments are evaluated more than once. o Customization When malloc() returns NULL, MAKE_1ARRAY prints a message on stderr and the program continues. The user can alter this behavior by redefining the MAKE_1ARRAY macro. For instance, to call exit() whenever malloc() fails, the user can do: #undef MAKE_1ARRAY #define MAKE_1ARRAY(a,n) do { \ (a) = malloc((n) * sizeof *(a)); \ if ((a)==NULL) { \ fprintf(stderr, "*** in file %s, function %s(), line %d: " \ "out of memory!\n", __FILE__, __func__, __LINE__); \ exit(EXIT_FAILURE); \ } \ } while (0) Since only MAKE_1ARRAY calls malloc() explicitly, this change affects the behavior of not only MAKE_1ARRAY but all other MAKE_* macros as well. ---- SAMPLE USAGE ------------- #include "array.h" int main(void) { float ***a; // can use any other type instead of "float" size_t p=3, q=4, r=5; // will make a 3x4x5 3-D array of float MAKE_3ARRAY(a, p, q, r); if (a==NULL) return EXIT_FAILURE; a[2][0][1] = 3.14; printf("%g \n", a[2][0][1]); FREE_3ARRAY(a); return EXIT_SUCCESS; } ---- END OF SAMPLE USAGE ------- Author: Rouben Rostamian, June 2000 Revised Jul 2000 Revised Sep 2002 Revised Jun 2003 - macros don't call exit anymore changed macro names to all-caps Revised Aug 2003 - changed reserved names to all caps Revised Dec 2003 - changed array index types to size_t (were int) - added an "out of memory" message, printed to stderr - most FREE* macros now come before MAKE* macros, for possible improved efficiency in preprocessing */ #ifndef H_ARRAY_ #define H_ARRAY_ #include #include /* ---------- 1D arrays ---------------------- */ #define MAKE_1ARRAY(a,n) do { \ (a) = malloc((n) * sizeof *(a)); \ if ((a)==NULL) \ fprintf(stderr, "*** in file %s, function %s(), line %d: " \ "out of memory!\n", __FILE__, __func__, __LINE__); \ } while (0) #define FREE_1ARRAY(a) do { \ free(a); \ a = NULL; \ } while (0) /* ---------- 2D arrays ---------------------- */ #define FREE_2ARRAY(a) do { \ size_t ARRAY_H2RESERVED; \ if (a==NULL) \ break; \ for (ARRAY_H2RESERVED=0; (a)[ARRAY_H2RESERVED]!=NULL; ARRAY_H2RESERVED++)\ FREE_1ARRAY((a)[ARRAY_H2RESERVED]); \ FREE_1ARRAY(a); \ } while (0) /* note: parenthesize first arg because it may be given as `*a' */ #define MAKE_2ARRAY(a,m,n) do { \ size_t ARRAY_H2RESERVED; \ MAKE_1ARRAY(a,(m)+1); \ if (a==NULL) \ break; \ (a)[m] = NULL; \ for (ARRAY_H2RESERVED=0; ARRAY_H2RESERVED<(m); ARRAY_H2RESERVED++) { \ MAKE_1ARRAY((a)[ARRAY_H2RESERVED],(n)); \ if ((a)[ARRAY_H2RESERVED]==NULL) { \ FREE_2ARRAY(a); \ break; \ } \ } \ } while (0) /* ---------- 3D arrays ---------------------- */ #define FREE_3ARRAY(a) do { \ size_t ARRAY_H3RESERVED; \ if (a==NULL) \ break; \ for (ARRAY_H3RESERVED=0; (a)[ARRAY_H3RESERVED]!=NULL; ARRAY_H3RESERVED++)\ FREE_2ARRAY((a)[ARRAY_H3RESERVED]); \ FREE_1ARRAY(a); \ } while (0) #define MAKE_3ARRAY(a,p,q,r) do { \ size_t ARRAY_H3RESERVED; \ MAKE_1ARRAY(a,(p)+1); \ if (a==NULL) \ break; \ (a)[p] = NULL; \ for (ARRAY_H3RESERVED=0; ARRAY_H3RESERVED<(p); ARRAY_H3RESERVED++) { \ MAKE_2ARRAY((a)[ARRAY_H3RESERVED],(q),(r)); \ if ((a)[ARRAY_H3RESERVED]==NULL) { \ FREE_3ARRAY(a); \ break; \ } \ } \ } while (0) /* ---------- 3D arrays ---------------------- */ #define FREE_4ARRAY(a) do { \ size_t ARRAY_H4RESERVED; \ if (a==NULL) \ break; \ for (ARRAY_H4RESERVED=0; (a)[ARRAY_H4RESERVED]!=NULL; ARRAY_H4RESERVED++)\ FREE_3ARRAY((a)[ARRAY_H4RESERVED]); \ FREE_1ARRAY(a); \ } while (0) #define MAKE_4ARRAY(a,p,q,r,s) do { \ size_t ARRAY_H4RESERVED; \ MAKE_1ARRAY(a,(p)+1); \ if (a==NULL) \ break; \ (a)[p] = NULL; \ for (ARRAY_H4RESERVED=0; ARRAY_H4RESERVED<(p); ARRAY_H4RESERVED++) { \ MAKE_3ARRAY((a)[ARRAY_H4RESERVED],(q),(r),(s)); \ if ((a)[ARRAY_H4RESERVED]==NULL) { \ FREE_4ARRAY(a); \ break; \ } \ } \ } while (0) /* ---------- synonyms ---------------------- */ #define MAKE_VECTOR(a,n) MAKE_1ARRAY(a,n) #define MAKE_MATRIX(a,m,n) MAKE_2ARRAY(a,m,n) #define FREE_VECTOR(a) FREE_1ARRAY(a) #define FREE_MATRIX(a) FREE_2ARRAY(a) #endif /* H_ARRAY_ */ CARP-1.0/QF.c0000664000076400007640000003047211356775462012523 0ustar volmelnvolmeln /* Author: Robert Davies, */ /* CARP: Clustering Algorithms' Referee Package Copyright (C) 2010 Volodymyr Melnykov and Ranjan Maitra This program is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version. This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details. You should have received a copy of the GNU General Public License along with this program. If not, see . Authors' contact information: Volodymyr Melnykov Ranjan Maitra Volodymyr.Melnykov@ndsu.edu maitra@iastate.edu Department of Statistics Department of Statistics North Dakota State University Iowa State University Fargo, ND 58102 Ames, IA 50011 */ #define UseDouble 0 /* all floating point double */ #include #include #include #include #define TRUE 1 #define FALSE 0 typedef int BOOL; #define pi 3.14159265358979 #define log28 .0866 /* log(2.0) / 8.0 */ static double sigsq, lmax, lmin, mean, c; static double intl, ersm; static int count, r, lim; static BOOL ndtsrt, fail; static int *n,*th; static double *lb,*nc; static jmp_buf env; static double exp1(double x) /* to avoid underflows */ { return x < -50.0 ? 0.0 : exp(x); } static void counter(void) /* count number of calls to errbd, truncation, cfe */ { count = count + 1; if ( count > lim ) longjmp(env,1); } static double square(double x) { return x*x; } static double cube(double x) { return x*x*x; } static double log1(double x, BOOL first) /* if (first) log(1 + x) ; else log(1 + x) - x */ { if (fabs(x) > 0.1) { return (first ? log(1.0 + x) : (log(1.0 + x) - x)); } else { double s, s1, term, y, k; y = x / (2.0 + x); term = 2.0 * cube(y); k = 3.0; s = (first ? 2.0 : - x) * y; y = square(y); for (s1 = s + term / k; s1 != s; s1 = s + term / k) { k = k + 2.0; term = term * y; s = s1; } return s; } } static void order(void) /* find order of absolute values of lb */ { int j, k; double lj; for ( j=0; j=0; k--) { if ( lj > fabs(lb[th[k]]) ) th[k + 1] = th[k]; else goto l1; } k = -1; l1 : th[k + 1] = j; } ndtsrt = FALSE; } static double errbd(double u, double* cx) /* find bound on tail probability using mgf, cutoff point returned to *cx */ { double sum1, lj, ncj, x, y, xconst; int j, nj; counter(); xconst = u * sigsq; sum1 = u * xconst; u = 2.0 * u; for (j=r-1; j>=0; j--) { nj = n[j]; lj = lb[j]; ncj = nc[j]; x = u * lj; y = 1.0 - x; xconst = xconst + lj * (ncj / y + nj) / y; sum1 = sum1 + ncj * square(x / y) + nj * (square(x) / y + log1(-x, FALSE )); } *cx = xconst; return exp1(-0.5 * sum1); } static double ctff(double accx, double* upn) /* find ctff so that p(qf > ctff) < accx if (upn > 0, p(qf < ctff) < accx otherwise */ { double u1, u2, u, rb, xconst, c1, c2; u2 = *upn; u1 = 0.0; c1 = mean; rb = 2.0 * ((u2 > 0.0) ? lmax : lmin); for (u = u2 / (1.0 + u2 * rb); errbd(u, &c2) > accx; u = u2 / (1.0 + u2 * rb)) { u1 = u2; c1 = c2; u2 = 2.0 * u2; } for (u = (c1 - mean) / (c2 - mean); u < 0.9; u = (c1 - mean) / (c2 - mean)) { u = (u1 + u2) / 2.0; if (errbd(u / (1.0 + u * rb), &xconst) > accx) { u1 = u; c1 = xconst; } else { u2 = u; c2 = xconst; } } *upn = u2; return c2; } static double truncation(double u, double tausq) /* bound integration error due to truncation at u */ { double sum1, sum2, prod1, prod2, prod3, lj, ncj, x, y, err1, err2; int j, nj, s; counter(); sum1 = 0.0; prod2 = 0.0; prod3 = 0.0; s = 0; sum2 = (sigsq + tausq) * square(u); prod1 = 2.0 * sum2; u = 2.0 * u; for (j=0; j 1.0) { prod2 = prod2 + nj * log(x); prod3 = prod3 + nj * log1(x, TRUE ); s = s + nj; } else prod1 = prod1 + nj * log1(x, TRUE ); } sum1 = 0.5 * sum1; prod2 = prod1 + prod2; prod3 = prod1 + prod3; x = exp1(-sum1 - 0.25 * prod2) / pi; y = exp1(-sum1 - 0.25 * prod3) / pi; err1 = ( s == 0 ) ? 1.0 : x * 2.0 / s; err2 = ( prod3 > 1.0 ) ? 2.5 * y : 1.0; if (err2 < err1) err1 = err2; x = 0.5 * sum2; err2 = ( x <= y ) ? 1.0 : y / x; return ( err1 < err2 ) ? err1 : err2; } static void findu(double* utx, double accx) /* find u such that truncation(u) < accx and truncation(u / 1.2) > accx */ { double u, ut; int i; static double divis[]={2.0,1.4,1.2,1.1}; ut = *utx; u = ut / 4.0; if ( truncation(u, 0.0) > accx ) { for ( u = ut; truncation(u, 0.0) > accx; u = ut) ut = ut * 4.0; } else { ut = u; for ( u = u / 4.0; truncation(u, 0.0) <= accx; u = u / 4.0 ) ut = u; } for ( i=0;i<4;i++) { u = ut/divis[i]; if ( truncation(u, 0.0) <= accx ) ut = u; } *utx = ut; } static void integrate(int nterm, double interv, double tausq, BOOL mainx) /* carry out integration with nterm terms, at stepsize interv. if (! mainx) multiply integrand by 1.0-exp(-0.5*tausq*u^2) */ { double inpi, u, sum1, sum2, sum3, x, y, z; int k, j, nj; inpi = interv / pi; for ( k = nterm; k>=0; k--) { u = (k + 0.5) * interv; sum1 = - 2.0 * u * c; sum2 = fabs(sum1); sum3 = - 0.5 * sigsq * square(u); for ( j = r-1; j>=0; j--) { nj = n[j]; x = 2.0 * lb[j] * u; y = square(x); sum3 = sum3 - 0.25 * nj * log1(y, TRUE ); y = nc[j] * x / (1.0 + y); z = nj * atan(x) + y; sum1 = sum1 + z; sum2 = sum2 + fabs(z); sum3 = sum3 - 0.5 * x * y; } x = inpi * exp1(sum3) / u; if ( ! mainx ) x = x * (1.0 - exp1(-0.5 * tausq * square(u))); sum1 = sin(0.5 * sum1) * x; sum2 = 0.5 * sum2 * x; intl = intl + sum1; ersm = ersm + sum2; } } static double cfe(double x) /* coef of tausq in error when convergence factor of exp1(-0.5*tausq*u^2) is used when df is evaluated at x */ { double axl, axl1, axl2, sxl, sum1, lj; int j, k, t; counter(); if (ndtsrt) order(); axl = fabs(x); sxl = (x>0.0) ? 1.0 : -1.0; sum1 = 0.0; for ( j = r-1; j>=0; j-- ) { t = th[j]; if ( lb[t] * sxl > 0.0 ) { lj = fabs(lb[t]); axl1 = axl - lj * (n[t] + nc[t]); axl2 = lj / log28; if ( axl1 > axl2 ) axl = axl1 ; else { if ( axl > axl2 ) axl = axl2; sum1 = (axl - axl1) / lj; for ( k = j-1; k>=0; k--) sum1 = sum1 + (n[th[k]] + nc[th[k]]); goto l; } } } l: if (sum1 > 100.0) { fail = TRUE; return 1.0; } else return pow(2.0,(sum1 / 4.0)) / (pi * square(axl)); } double qfc(double* lb1, double* nc1, int* n1, int *r1in, double *sigmain, double *c1in, int *lim1in, double *accin, double* trace, int* ifault) /* distribution function of a linear combination of non-central chi-squared random variables : input: lb[j] coefficient of j-th chi-squared variable nc[j] non-centrality parameter n[j] degrees of freedom j = 0, 2 ... r-1 sigma coefficient of standard normal variable c point at which df is to be evaluated lim maximum number of terms in integration acc maximum error output: ifault = 1 required accuracy NOT achieved 2 round-off error possibly significant 3 invalid parameters 4 unable to locate integration parameters 5 out of memory trace[0] absolute sum trace[1] total number of integration terms trace[2] number of integrations trace[3] integration interval in final integration trace[4] truncation point in initial integration trace[5] s.d. of initial convergence factor trace[6] cycles to locate integration parameters */ { int r1 = r1in[0], lim1 = lim1in[0]; double sigma = sigmain[0], c1 = c1in[0], acc = accin[0]; int j, nj, nt, ntm; double acc1, almx, xlim, xnt, xntm; double utx, tausq, sd, intv, intv1, x, up, un, d1, d2, lj, ncj; double qfval = 0; static int rats[]={1,2,4,8}; if (setjmp(env) != 0) { *ifault=4; goto endofproc; } r=r1; lim=lim1; c=c1; n=n1; lb=lb1; nc=nc1; for ( j = 0; j<7; j++ ) trace[j] = 0.0; *ifault = 0; count = 0; intl = 0.0; ersm = 0.0; qfval = -1.0; acc1 = acc; ndtsrt = TRUE; fail = FALSE; xlim = (double)lim; th=(int*)malloc(r*(sizeof(int))); if (! th) { *ifault=5; goto endofproc; } /* find mean, sd, max and min of lb, check that parameter values are valid */ sigsq = square(sigma); sd = sigsq; lmax = 0.0; lmin = 0.0; mean = 0.0; for (j=0; j lj) lmin = lj; } if ( sd == 0.0 ) { qfval = (c > 0.0) ? 1.0 : 0.0; goto endofproc; } if ( lmin == 0.0 && lmax == 0.0 && sigma == 0.0 ) { *ifault = 3; goto endofproc; } sd = sqrt(sd); almx = (lmax < - lmin) ? - lmin : lmax; /* starting values for findu, ctff */ utx = 16.0 / sd; up = 4.5 / sd; un = - up; /* truncation point with no convergence factor */ findu(&utx, .5 * acc1); /* does convergence factor help */ if (c != 0.0 && (almx > 0.07 * sd)) { tausq = .25 * acc1 / cfe(c); if (fail) fail = FALSE ; else if (truncation(utx, tausq) < .2 * acc1) { sigsq = sigsq + tausq; findu(&utx, .25 * acc1); trace[5] = sqrt(tausq); } } trace[4] = utx; acc1 = 0.5 * acc1; /* find RANGE of distribution, quit if outside this */ l1: d1 = ctff(acc1, &up) - c; if (d1 < 0.0) { qfval = 1.0; goto endofproc; } d2 = c - ctff(acc1, &un); if (d2 < 0.0) { qfval = 0.0; goto endofproc; } /* find integration interval */ intv = 2.0 * pi / ((d1 > d2) ? d1 : d2); /* calculate number of terms required for main and auxillary integrations */ xnt = utx / intv; xntm = 3.0 / sqrt(acc1); if (xnt > xntm * 1.5) { /* parameters for auxillary integration */ if (xntm > xlim) { *ifault = 1; goto endofproc; } ntm = (int)floor(xntm+0.5); intv1 = utx / ntm; x = 2.0 * pi / intv1; if (x <= fabs(c)) goto l2; /* calculate convergence factor */ tausq = .33 * acc1 / (1.1 * (cfe(c - x) + cfe(c + x))); if (fail) goto l2; acc1 = .67 * acc1; /* auxillary integration */ integrate(ntm, intv1, tausq, FALSE ); xlim = xlim - xntm; sigsq = sigsq + tausq; trace[2] = trace[2] + 1; trace[1] = trace[1] + ntm + 1; /* find truncation point with new convergence factor */ findu(&utx, .25 * acc1); acc1 = 0.75 * acc1; goto l1; } /* main integration */ l2: trace[3] = intv; if (xnt > xlim) { *ifault = 1; goto endofproc; } nt = (int)floor(xnt+0.5); integrate(nt, intv, 0.0, TRUE ); trace[2] = trace[2] + 1; trace[1] = trace[1] + nt + 1; qfval = 0.5 - intl; trace[0] = ersm; /* test whether round-off error could be significant allow for radix 8 or 16 machines */ up=ersm; x = up + acc / 10.0; for (j=0;j<4;j++) { if (rats[j] * x == rats[j] * up) *ifault = 2; } endofproc : free((char*)th); trace[6] = (double)count; return qfval; } CARP-1.0/hierclust.c0000664000076400007640000002566711357661203014215 0ustar volmelnvolmeln /* CARP: Clustering Algorithms' Referee Package Copyright (C) 2010 Volodymyr Melnykov and Ranjan Maitra This program is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version. This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details. You should have received a copy of the GNU General Public License along with this program. If not, see . Authors' contact information: Volodymyr Melnykov Ranjan Maitra Volodymyr.Melnykov@ndsu.edu maitra@iastate.edu Department of Statistics Department of Statistics North Dakota State University Iowa State University Fargo, ND 58102 Ames, IA 50011 */ #include #include #include #include #include "array.h" #define inf 1e+40; /* Given a hierarchical clustering through a sequence of agglomerations, derive the assignment for the top (lev-1)th level of the hierarchy. Parameters: n: number of observations ia,ib: vectors of dimension N defining the agglomerations lev: number of clusters in partition. iclass: n-dimensional vector of cluster assignments C code written by Ranjan Maitra, Baltimore (07/02/02) */ void hclass(int n, int *ia, int *ib, int lev, int *iclass) { int i,j,k; for (i=0;i=0;k--) { if (iclass[ia[k]]==j) { iclass[ib[k]]=j; } } j--; } j=0; iclass[ia[n-2]]=j; for (k=(n-lev-1);k>=0;k--) { if (iclass[ia[k]]==j) { iclass[ib[k]]=j; } } return; } /* HIERARCHICAL CLUSTERING using (user-specified) criterion. Parameters: n,m number of observations and coordinates data(n,m) input data matrix, iopt clustering criterion to be used, IA, IB, CRIT history of agglomerations; dimension N-1 each C code rewritten by Ranjan Maitra, Baltimore (07/02/02) from f77 code of F. Murtagh, ESA/ESO/STECF, Garching (1986). Note that I use the true distances, not squared Euclidean which FM used. */ int ioffset(int n, int i, int j) { /* Map (i,j)th element of upper half */ return(j+i*n-(i+1)*(i+2)/2); /* diagonal symmetric matrix to vector */ } int imin(int aa, int bb) { if (aa1) { /* Next, determine least diss. using list of NNs */ int jj=0,jm=0,im=0; dmin=inf; for (i= 0;i<(n-1);++i) { if (flag[i]==1) { if (disnn[i]. Authors' contact information: Volodymyr Melnykov Ranjan Maitra Volodymyr.Melnykov@ndsu.edu maitra@iastate.edu Department of Statistics Department of Statistics North Dakota State University Iowa State University Fargo, ND 58102 Ames, IA 50011 */ /* This is the main program that reads parameters and runs the following 3 steps: 1. C-MixSim which generates mixtures with prespecified level of overlap and datasets from these mixtures. 2. Clustering algorithm provided by user. 3. AdjRand which computes Adjusted Rand index for simulated datasets. */ #include #include #include #include #include #include int getopt(int argc, char *const argv[], const char *optstring); extern char *optarg; extern int optind, optopt; int main(int argc, char *argv[]) { if (argc == 1){ printf("Usage:\n> ./CARP \n\nParameters:\n"); printf("-b: average overlap (no default value)\n"); printf("-m: maximum overlap (no default value)\n"); printf("-p: number of dimensions (2 by default)\n"); printf("-K : number of mixing components (2 by default)\n"); printf("-n: number of observations generated from every mixture (0 by default)\n"); printf("-#: number of simulated mixtures (1 by default)\n"); printf("-0 : name of clustering program (no default name)\n"); printf("-1 : name of partition analyzing program (“AdjRand” by default)\n"); printf("-s: spherical covariance matrix structure (non-spherical by default if option is unspecified)\n"); printf("-e: maximum eccentricity (0.90 by default)\n"); printf("-z : smallest mixing proportion (equal mixing proportions 1/K by default)\n"); printf("-u: upper bound for Uniform(0, upper-bound ) distribution from which mean vectors are generated\n"); printf("-r : maximum number of resimulations (100 by default)\n"); printf("-a: accuracy of estimation (1e-06 by default)\n"); printf("-l : maximum number of integration terms (1e06 by default)\n"); printf("-P : name of the file containing mixing proportions (“Pi.dat” by default)\n"); printf("-M : name of the file containing mean vectors (“Mu.dat” by default)\n"); printf("-S : name of the file containing covariance matrices in triangular form (“LTSigma.dat” by default)\n"); printf("-D: name of working directory (“DATA” by default)\n"); printf("-I : name of the file containing numbers of observations generated from every cluster (“Nk.dat” by default)\n"); printf("-i : name of the file containing estimated classifications (“idEst.dat” by default)\n"); printf("-X : name of the file containing simulated datasets (“x.dat” by default)\n"); printf("-W : name of the file containing maps of pairwise overlaps (“overMap.dat” by default)\n"); printf("-C : name of the file containing characteristics of simulated mixtures (“overBarMax.dat” by default)\n"); printf("-R: name of the file containing index values (Adjusted Rand index and “AR.dat” by default)\n\n"); exit(1); } int c, errflag0, errflag1, code; char file1[500] = "./C-MixSim"; char file2[500] = "./"; char file3[500] = "./"; char options1[500] = ""; char options2[500] = ""; char options3[500] = ""; errflag0 = 1; errflag1 = 1; /* DEALING WITH PARAMETERS */ while ((c = getopt(argc, argv, ":b:m:p:K:se:z:u:r:a:l:P:M:S:D:I:i:X:W:C:R:n:#:0:1:")) != -1) { switch(c) { case '0': strcat(file2, optarg); errflag0 = 0; break; case '1': strcat(file3, optarg); errflag1 = 0; break; default: if ((c != '0') & (c != 'i') & (c != 'R')){ strcat(options1, " "); strcat(options1, argv[optind-1]); } if ((c == 'p') | (c == 'K') | (c == 'n') | (c == '#') | (c == 'D') | (c == 'X') | (c == 'i')){ strcat(options2, " "); strcat(options2, argv[optind-1]); } if ((c == 'K') | (c == 'n') | (c == '#') | (c == 'D') | (c == 'I') | (c == 'i') | (c == 'R')){ strcat(options3, " "); strcat(options3, argv[optind-1]); } break; } } /* default analyzing program */ if (errflag1 != 0){ strcat(file3, "AdjRand"); } strcat(file1, options1); strcat(file2, options2); strcat(file3, options3); // if (errflag0 != 0){ // printf("ERROR: No clustering program has been detected [option -0]\n"); // exit(1); // } /* STEP 1 - Mixture simulation */ printf("\n\nSTEP 1: Running %s\n\n", file1); code = system(file1); if (code != 0){ printf("ERROR: Unable to fulfill STEP 1\n"); exit(1); } if (errflag0 == 0){ /* STEP 2 - Run of clustering algorithm */ printf("\n\n\nSTEP 2: Running %s\n\n\n", file2); code = system(file2); if (code != 0){ printf("ERROR: Unable to fulfill STEP 2\n"); exit(1); } /* STEP 3 - Calculation of Adjusted Rand index */ printf("\n\n\nSTEP 3: Running %s\n\n\n", file3); code = system(file3); if (code != 0){ printf("ERROR: Unable to fulfill STEP 3\n"); exit(1); } } return 0; } CARP-1.0/minunit.h0000644000076400007640000000236111356774723013676 0ustar volmelnvolmeln /* CARP: Clustering Algorithms' Referee Package Copyright (C) 2010 Volodymyr Melnykov and Ranjan Maitra This program is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version. This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details. You should have received a copy of the GNU General Public License along with this program. If not, see . Authors' contact information: Volodymyr Melnykov Ranjan Maitra Volodymyr.Melnykov@ndsu.edu maitra@iastate.edu Department of Statistics Department of Statistics North Dakota State University Iowa State University Fargo, ND 58102 Ames, IA 50011 */ /* file: minunit.h */ #ifndef MINUNIT_H #define MINUNIT_H #define func_test(test, message1, message2) do { \ if (test) printf(message1); \ else printf(message2); \ } while (0) #endif /* MINUNIT_H */ CARP-1.0/TEST/0000755000076400007640000000000011357676030012610 5ustar volmelnvolmelnCARP-1.0/TEST/Mu.test0000664000076400007640000000071611324705272014073 0ustar volmelnvolmeln1.252247e-03 1.210485e-01 9.369321e-01 6.814153e-01 3.689626e-01 1.587390e-01 6.093170e-01 9.835994e-02 5.947354e-01 1.960200e-01 3.310922e-01 8.731364e-01 3.875023e-01 9.396728e-01 1.252729e-01 1.106599e-01 3.865950e-01 2.148268e-01 9.270482e-01 4.715763e-02 1.921139e-01 6.903454e-01 6.311570e-01 5.949675e-01 7.657742e-01 3.360492e-01 8.787721e-01 1.807687e-02 9.408491e-02 5.964270e-02 1.882026e-01 9.533716e-02 1.806912e-01 1.251347e-01 7.767524e-01 CARP-1.0/TEST/Pi.test0000664000076400007640000000013411353171471014054 0ustar volmelnvolmeln1.831131e-01 1.112458e-01 1.405014e-01 1.442849e-01 1.089084e-01 1.549885e-01 1.569579e-01 CARP-1.0/TEST/Pi.iris0000666000076400007640000000005211357473034014050 0ustar volmelnvolmeln3.3333333e-01 3.3333333e-01 3.3333333e-01 CARP-1.0/TEST/Mu.iris0000666000076400007640000000011211357471157014062 0ustar volmelnvolmeln5.006 3.428 1.462 0.246 5.936 2.770 4.260 1.326 6.588 2.974 5.552 2.026 CARP-1.0/TEST/.RData0000755000076400007640000003406111263754702013612 0ustar volmelnvolmeln|7ۿ$"hFZD!KE%*"QiQR*UIlmiRe`m>f̠"-~>999>'^u^u/h8x`999y9y%ӟ}v KVe\yֵe.NOݫŲM }GgN _2-58~jm kΑeL+uI5~nMSueE&R`1L^3>X{(pDoSprGs'ضfA&ϝ ߥNSF=%䄼3C΅unB㌝)'mlI͡kFCuT\sx85c2\cE`yض{j^l pZ7FxgyvSx+6Up<r+ahOg?3 ++ܢo3dAIx`{Y_>Hf p?{,q1c}knm.شjk) ۟^|KΏ;hejDhh8Ku8Gp>4Ƃ\8S%|7erӿg ''c_%˂1`M'gc0xU ԋG̘U7*o v=ux SֱN<"Q23s4_{w@p&yN+ MzC3ĎgI* f `nwwi68o&ιFl`wtޅNgRRxD"5'##iS 7? C?7~? pzy+ ' }`.a#nY:w;;{tװ{<.Fѝ=Yg*<~lӎsMʲmo=tO-$6~b.| 瓖8&L"ύ;i̝c[*ޑ[Eu1Y̊[h SR =A[7|ne9ֲ,EGs0 p$z/=lOAw V=vŜO*>Kx!025K4ǃg"}0oΖ$;g?o<\&ws BׁǸ;bA0UZxjy㶹jo?9y$]p_w+u7|/m#yL~S6%sB?~~ROxS35rW@8rp,h, as֭@O> +|+=~ /vXۭsx#.EG vd=]]8M;?C4L;znuWh,cf!"A.E: Xv&eV=M?N@lo7߯#qKA`]t \sIxT/ R;u:r1~XYO΋XSfROǠ`sp${s\NN_S\\94s>s{%XgdGBhHOfS Oo[ B]sX x6}Q-Y~> ۂ;$gVh >L1ML ܿ!rdP>xuhmա;zFdpNK(oo{8bYp~|$'lN:9 )x+2׃Cݼ `#fA]ē-βbr mkkշڮYQ! 8_/ " f>oWĀ'^k}ќ-$&h*]6|J3躄5?il\ 5$'L(h[m}qKzi?R2ۦkL<zB[pqy{&%YqqO_2Hڜo"_p7ዧR~YOepV^Y>Tε UWUiZ5Ε#6qqt\.{w<}@fSz ؍6c^}0cfٽɛWA+2+ mrW A¶ v=6]lj\?!: O/%mYaU_:zӛfH|l8xEIm41w?r~^45Yz^':PXT>𵏹>t:ndI(bmPN%/ %P3_J\pԷ~6\LKJн]tFj,z%fK5_46UK 5oL1.l(_$,'T\>N]b1oQvt[U  έF!Cʁ jMϠ27S.~(T\v]Ĉe QƻGO(VOqKcPg=gow!*2l}{ J$zͷQ2|7f&(a-e֣:O5n$j2 ۰ߗxstLl^~cl8ͦwY1F?hLYuFvtF63N!& `m!^`6 H54X-zPA&rɊ:ӨM.(PYߡ~*6rGQ *gzFu;'*Dk,Sj_Fw"`jTxKz ,^1c=bQo q/ԏT&r̠'U'P"2?_?6-Gm uVhA岊Tem/zlh8R9~zJK2ϹnG6{/;z;wN^I(!zyhݧew t.LqyjNbBSQ#o8V֫w9] b+Q7)hs*g,Լ |(G%SZKd30N:EOGdWӲ8Tę{\C9%٨ҝ4g͑y ,T* *]KޓgjFmEQA[i Z [P9Q`^ wi+&tp^:6XN_MW T)h]Q 5g|ڱ̰(Ч- Dڍ:Ơf/Q/7Wի](.CoA!f$>杪ˢP58co 堊cTjYz$ϛ(*C%% ǽKcI-Ъ/QIyjE>V_ ?{n %؁tΥwzDo]t D7w'M/ ^˂ȮƊʍ7]㭩dO 訙BkO|"ܙan-^,'²wo眅h ^}U} Yă Xp3z"O!,gDC5o?]O0'Q« '@c!y3Шj^&C 4o:+<-дXvY GSdPyNѕ{&.#$`Si$]>ؔ]Y]3YD`Z* j}mmc;4M޲:B/s?'47!2-msv5 [Aս!Mvx_{$ǀ|ՁQ^ʁPqR\ $ק/":,!˞W YLE9}`ts'306%RfnXk! BK^D~Òoݶg@Z"(?oƟ͎ c4Cpn`ӅBS4O˃rM;%(>QB1ѵvip҆`ͳ23e W@$_ă3 s_N>D x}i:Dޯڿ?i|""q(,eBG93,YaOXaO9lVAK_~͉}k :8 wiU!6UlM`Q yC SE2j_z!VY}$6(hC$Peِ3BcѼ[9/e{ $`j4 {B|J{u7XBRNAc>hxa*$3?I0˾ e8 yݾ o~qw[hӭMg{bǮ k+~(@؟"C7Ҽ˽{c!~bH˚ZqmF4xU ͈MS!I@Ajk UvCvɐnɌxIeU/_Lw#^wAr~^ѹG5N hQ3(ĠnjckO= ECoM &hhHV^_ 3$JT=F+@ܶ%ؓb޼HC|fwHr)CЍ|ά=qj yAў2 aOCrՒ8h=LJP [ѼU+]I{;Uy61oix$ĿG Hk(& IVmrVL1ğ ;邦p cD{]\ZBJ9OlX>I7]8fb ozߐivX(߹ʄwͫKoy͊h]C m ~,E0*݌G}~ڛQza/QFiQ>ECb׮dĀVg__򹾼F<תT,: qFo߂$yiL.BV~+wMwYe4uxZ 7^ʼnʅ׎ĴUmd;778ރ)ROY›YdmA2j*jIֹ{nzNm.ڗl$Qhs\o %ӣ?@"1[!)3 H>re_t5xl9t⸁r@:M]uN6&8&BS*E}^r 9ZZGe/I9yd~mt "~5_UFSM( I€~<ⴅF 8>2D>Q7bc@^!mL!r񮳷z&9.oĻu؎&f?^C+tΊ6>ߝOZAڗ.w'u$cBYӓyo~/$^%(?(ku?x#Zd?ϛxاVH"ơ`H&,Gnѯ'؎:(,NFpn,q4`g2ĢF@;[?yi=pU]dFQ;kB,5%ڞc[RmHz:w&չ ~JMɹ9Jw >K3ĦDkA'CI /4w]| ᱐W6ho#h۽+ȋK{U jNp,΁6KT2&\C=N^7W8G gHmuxֶ`shN/?+CR;Y+Tr:>nR9)/@l*i@ Rktss_+FK] i&_ yH%X ~!; %Y7MU;7 niXdȢps6hTy=o\w ^Z9*CaBAaTZt\sBTtcdԩBED7dɕ۳Bfmv2-Ct ;GBV4_%% Y2ʏH0ևL{Ŝq6],;h|ɐw btHhNh5G; }Pt+ǣEw]K@6*#!+9Q/d@JA\*D]jA >L6 }Kh>d=\]S!UHC9xw8ђķGI|:|Y9 'U ˃^E Օf'{7#./BV yne`/_ ՜>|~`GbͰğ7|zZ>NJXo^فy{HfR@ d:2uA꤁JUp %9 s&BܫdrP-QE/BnIm'ٷT.d]뼈.sHS'd*/L.Vq˻:׈ʇyw-.?;EDCWeZZ{ hc3<> }H|Zx⟹we!EoZO*덫2Ur<*h&v?sChὠ: c[WoTrU"rODLEsodCvLD%ZzHt Y-/1FDEh_v)H+rB08d5%qtX'@66<ï3/}ܬ{{i]r{"({ɶ-ItB4 3#kv"GBmˈL<5&C,H//X>2 >7B}m=ZҗW>sӮ-_!>Hpgo?T,Qf2-/gFK_LFN0&_6 "Щtrŀ*|L_z@:oh}m7n9Asa([hQ_:!Z[ւK8#H8+C.Șf -Q`4CŧWyG*oGwRP4^7֓fuD6b5fUE"F#Fm3}J Sk1abP6=ԗOBwP6|EߗrF6nƑ(ة_@ p $)^|i4wwyý\Fcy(VQ(F){xh磬2: QP}[(UW59# #\>(z2$庫x(b<r;QBI;}V q椛i( 1Z|Fr ݠ^7Q8JPR)BAq^P*.OaX>G;U@ )ysK亸YͶ љ:,~ȧJ@#j%8ltXoU>hai8g<`dƭ6sk"G=k (mRށٍٶ(i盉C6Lg.Bم}Q@\?d6l9Co^$oM7QׯRLӬm9Vxy=U(^adg1E_{5E޵RL2نsZ^3Z~k$F~=t$YH9x}L+s^ODA]RY,|1RƣhSC.^|Z?&I*8< /^0~ de^znOAIBqкoG7FЋIr/OPV~k%TJ0袃>R8-.|tm`,h3|xMG@!FU_wcE5iچIivGZ%EQNe#'8<@ٴ2EW(~DX9t]ܘ!ylSakB r֑?JdlEEГlPOJ;j4E~q"v]{GݥQڧWrZuvY-"E+&h sm>Fk7P0ҩ(N7Q@ywxD GΑQ;Au7lsdІ}UOr-rB]y6g/Ikkһ? m:E&m*QB JP@(lk.םNţ ;A[pFyP,O(tzaȩ(z*P21Y/(2`3KUbA;B+/*(ZWclMgMP؇C忨/TnK} 3OBE~ߧB3-3Bn y;z?NGsƓ:^QE}n?xJff; b-w0tȍ~{ i zw Je?Dr֙{W?gPMF hq/@yOY40ض=LϵBүbse_1ky8 *%W$dmUSޠj;QaZr' dZfo^#uPuaʩ$(DeQLC -fPr-(Pd_;W{}zyP8ޱz+:HަyLEyèpd=v>jx'JiĢL CB/*vw{g5k*#JN@pq@uY0`Q ?{b8;En3ΔG^b Ww%@e}%TMV+;U?5OFIvkK=W'*w^Fů[m{f-QpgXn2*m6o]t6Oa[T}/Qׇ]J%xIUo̞fcBskX<;0m6݌SFO ѵ7лRxAZ*ϕ>+e+n;J;̐1jS#e0. Ioūy]yDT&5ωQBPyN5%/~=;Q6&8,*?Bԛe% ]30Azc!*N] _E먭2z~,Ng5: y<vF.Q*՝EC@և2ݼ;(.ʿ^ͷ1ʵ6N 7kv@}߽l0gy-F.JPߡҷ:cPO7: fp+sO'yEPcUHDGu`J-PGNfoqG'无.PZ-Cmx'>x4_fө YXкw++%5:Qq6;tN)ƥ3'CEOzUU.1c̓L~HA{ uǻs?*us"@] }_ >TR(vfwT0dCToTgs_{}vK﹦ $0RX`FP sݺ3P53lM"~įn8jAqR޷Qw_F.[l(-9; ? B%rII@9L{uIҊơv YqoԍzJ(U[' D1ל~XqHTK70XPL]J3mU\љŨ8|>,fS>B~TaE@;%DpN]vڰTW/dSBY~C>Tl4b\ϖW~@5'> 9ʾ{0t`O=l#*=mT8 LݽFEOφ~yoY:%)5?m5Kݞx< .(&p?|OoT 3M U)AJ?"~O L !2^F]jJ}֘StA`9+Pux?}'Q\&W䧛^H^CARP-1.0/TEST/LTSigma.test0000664000076400007640000000257711324705272015021 0ustar volmelnvolmeln1.621167e-02 4.687798e-03 1.783813e-02 4.930906e-04 1.637381e-03 1.819005e-02 -6.684429e-03 -3.113320e-03 -4.424269e-03 1.340989e-02 1.040291e-02 -4.613393e-04 -4.847592e-03 -6.120238e-03 2.381067e-02 1.621573e-02 -5.858855e-04 3.362437e-02 8.637641e-04 1.071124e-02 1.770941e-02 8.434695e-04 5.302939e-03 6.407430e-03 3.707597e-02 -5.479690e-03 1.826564e-02 5.290868e-03 4.672782e-04 3.467262e-02 2.498838e-02 -1.092985e-03 3.074239e-02 -6.018394e-03 9.633310e-03 2.438968e-02 -1.505005e-03 1.300912e-02 7.436055e-03 1.959343e-02 4.244509e-03 -7.281346e-03 -1.164089e-02 -6.154002e-03 1.999118e-02 1.821750e-02 -8.124875e-03 3.350045e-02 3.134526e-03 -2.118332e-02 4.042105e-02 -3.561823e-03 -3.153373e-03 5.905472e-03 4.587827e-02 -4.717958e-03 -4.866358e-05 8.449999e-03 -1.850409e-03 2.193378e-02 2.568136e-02 -7.965512e-03 2.076638e-02 -7.996792e-03 1.746973e-03 1.134245e-02 4.745011e-03 -1.590397e-03 -3.289079e-03 1.037866e-02 -5.133851e-03 2.492457e-03 1.554633e-03 -4.796974e-04 9.017570e-03 1.201213e-02 -1.864540e-03 1.003660e-02 2.490482e-03 2.063010e-03 8.533450e-03 2.276079e-03 1.518940e-03 4.118088e-03 2.444800e-02 8.475646e-04 -2.864911e-03 -1.748466e-03 -5.619680e-03 8.362936e-03 6.472221e-03 2.552347e-04 1.157326e-02 1.256301e-03 5.321933e-03 1.750908e-02 1.191448e-03 -2.810004e-04 2.407939e-03 9.799328e-03 6.570816e-04 -3.569637e-03 -9.303376e-03 -3.481258e-04 2.087818e-02 CARP-1.0/TEST/.Rhistory0000755000076400007640000000057311263754702014443 0ustar volmelnvolmelnx = matrix(scan("x.dat"), ncol = 3, byrow = T) head(x) plot(x[,2], x[,3]) plot(x[,2], x[,3], cex = 0.5, pch = 19) x = matrix(scan("x.dat"), ncol = 3, byrow = T) plot(x[,2], x[,3], cex = 0.5, pch = 19) x = matrix(scan("x.dat"), ncol = 3, byrow = T) plot(x[,2], x[,3], cex = 0.5, pch = 19) x = matrix(scan("x.dat"), ncol = 3, byrow = T) plot(x[,2], x[,3], cex = 0.5, pch = 19) q() CARP-1.0/TEST/LTSigma.iris0000666000076400007640000000051711357475070015010 0ustar volmelnvolmeln0.12424898 0.09921633 0.14368980 0.01635510 0.01169796 0.030159184 0.01033061 0.00929796 0.006069388 0.011106122 0.26643265 0.08518367 0.09846939 0.18289796 0.08265306 0.22081633 0.05577959 0.04120408 0.07310204 0.03910612 0.40434286 0.09376327 0.10400408 0.30328980 0.07137959 0.30458776 0.04909388 0.04762857 0.04882449 0.07543265 CARP-1.0/cephes_eigens.h0000644000076400007640000000224011356774520015003 0ustar volmelnvolmeln /* CARP: Clustering Algorithms' Referee Package Copyright (C) 2010 Volodymyr Melnykov and Ranjan Maitra This program is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version. This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details. You should have received a copy of the GNU General Public License along with this program. If not, see . Authors' contact information: Volodymyr Melnykov Ranjan Maitra Volodymyr.Melnykov@ndsu.edu maitra@iastate.edu Department of Statistics Department of Statistics North Dakota State University Iowa State University Fargo, ND 58102 Ames, IA 50011 */ #ifndef CEPHES_EIGENS_H #define CEPHES_EIGENS_H void cephes_eigens(double *A, double *EV, double *E, int N); #endif /* CEPHES_EIGENS_H */ CARP-1.0/makefile0000664000076400007640000000074311356527253013540 0ustar volmelnvolmelnCC = gcc CFLAGS = -Wall -O3 -std=c99 -pedantic OBJ = cephes_eigens.o GenData.o GenPars.o Matrix.o Overlap.o Print.o QF.o Random.o symmEigValDec.o CARP: CARP.c $(OBJ) $(CC) $(CFLAGS) CARP.c -o CARP $(OBJ) -lm $(CC) $(CFLAGS) C-MixSim.c -o C-MixSim $(OBJ) -lm $(CC) $(CFLAGS) AdjRand.c -o AdjRand -lm $(CC) $(CFLAGS) hierclust.c -o hierclust -lm $(CC) $(CFLAGS) unittest.c -o unittest $(OBJ) -lm clean: rm -rf *.o CARP C-MixSim AdjRand hierclust unittest check: ./unittest CARP-1.0/DATA/0000755000076400007640000000000011357707117012543 5ustar volmelnvolmelnCARP-1.0/DATA/.RData0000755000076400007640000003406111263754702013544 0ustar volmelnvolmeln|7ۿ$"hFZD!KE%*"QiQR*UIlmiRe`m>f̠"-~>999>'^u^u/h8x`999y9y%ӟ}v KVe\yֵe.NOݫŲM }GgN _2-58~jm kΑeL+uI5~nMSueE&R`1L^3>X{(pDoSprGs'ضfA&ϝ ߥNSF=%䄼3C΅unB㌝)'mlI͡kFCuT\sx85c2\cE`yض{j^l pZ7FxgyvSx+6Up<r+ahOg?3 ++ܢo3dAIx`{Y_>Hf p?{,q1c}knm.شjk) ۟^|KΏ;hejDhh8Ku8Gp>4Ƃ\8S%|7erӿg ''c_%˂1`M'gc0xU ԋG̘U7*o v=ux SֱN<"Q23s4_{w@p&yN+ MzC3ĎgI* f `nwwi68o&ιFl`wtޅNgRRxD"5'##iS 7? C?7~? pzy+ ' }`.a#nY:w;;{tװ{<.Fѝ=Yg*<~lӎsMʲmo=tO-$6~b.| 瓖8&L"ύ;i̝c[*ޑ[Eu1Y̊[h SR =A[7|ne9ֲ,EGs0 p$z/=lOAw V=vŜO*>Kx!025K4ǃg"}0oΖ$;g?o<\&ws BׁǸ;bA0UZxjy㶹jo?9y$]p_w+u7|/m#yL~S6%sB?~~ROxS35rW@8rp,h, as֭@O> +|+=~ /vXۭsx#.EG vd=]]8M;?C4L;znuWh,cf!"A.E: Xv&eV=M?N@lo7߯#qKA`]t \sIxT/ R;u:r1~XYO΋XSfROǠ`sp${s\NN_S\\94s>s{%XgdGBhHOfS Oo[ B]sX x6}Q-Y~> ۂ;$gVh >L1ML ܿ!rdP>xuhmա;zFdpNK(oo{8bYp~|$'lN:9 )x+2׃Cݼ `#fA]ē-βbr mkkշڮYQ! 8_/ " f>oWĀ'^k}ќ-$&h*]6|J3躄5?il\ 5$'L(h[m}qKzi?R2ۦkL<zB[pqy{&%YqqO_2Hڜo"_p7ዧR~YOepV^Y>Tε UWUiZ5Ε#6qqt\.{w<}@fSz ؍6c^}0cfٽɛWA+2+ mrW A¶ v=6]lj\?!: O/%mYaU_:zӛfH|l8xEIm41w?r~^45Yz^':PXT>𵏹>t:ndI(bmPN%/ %P3_J\pԷ~6\LKJн]tFj,z%fK5_46UK 5oL1.l(_$,'T\>N]b1oQvt[U  έF!Cʁ jMϠ27S.~(T\v]Ĉe QƻGO(VOqKcPg=gow!*2l}{ J$zͷQ2|7f&(a-e֣:O5n$j2 ۰ߗxstLl^~cl8ͦwY1F?hLYuFvtF63N!& `m!^`6 H54X-zPA&rɊ:ӨM.(PYߡ~*6rGQ *gzFu;'*Dk,Sj_Fw"`jTxKz ,^1c=bQo q/ԏT&r̠'U'P"2?_?6-Gm uVhA岊Tem/zlh8R9~zJK2ϹnG6{/;z;wN^I(!zyhݧew t.LqyjNbBSQ#o8V֫w9] b+Q7)hs*g,Լ |(G%SZKd30N:EOGdWӲ8Tę{\C9%٨ҝ4g͑y ,T* *]KޓgjFmEQA[i Z [P9Q`^ wi+&tp^:6XN_MW T)h]Q 5g|ڱ̰(Ч- Dڍ:Ơf/Q/7Wի](.CoA!f$>杪ˢP58co 堊cTjYz$ϛ(*C%% ǽKcI-Ъ/QIyjE>V_ ?{n %؁tΥwzDo]t D7w'M/ ^˂ȮƊʍ7]㭩dO 訙BkO|"ܙan-^,'²wo眅h ^}U} Yă Xp3z"O!,gDC5o?]O0'Q« '@c!y3Шj^&C 4o:+<-дXvY GSdPyNѕ{&.#$`Si$]>ؔ]Y]3YD`Z* j}mmc;4M޲:B/s?'47!2-msv5 [Aս!Mvx_{$ǀ|ՁQ^ʁPqR\ $ק/":,!˞W YLE9}`ts'306%RfnXk! BK^D~Òoݶg@Z"(?oƟ͎ c4Cpn`ӅBS4O˃rM;%(>QB1ѵvip҆`ͳ23e W@$_ă3 s_N>D x}i:Dޯڿ?i|""q(,eBG93,YaOXaO9lVAK_~͉}k :8 wiU!6UlM`Q yC SE2j_z!VY}$6(hC$Peِ3BcѼ[9/e{ $`j4 {B|J{u7XBRNAc>hxa*$3?I0˾ e8 yݾ o~qw[hӭMg{bǮ k+~(@؟"C7Ҽ˽{c!~bH˚ZqmF4xU ͈MS!I@Ajk UvCvɐnɌxIeU/_Lw#^wAr~^ѹG5N hQ3(ĠnjckO= ECoM &hhHV^_ 3$JT=F+@ܶ%ؓb޼HC|fwHr)CЍ|ά=qj yAў2 aOCrՒ8h=LJP [ѼU+]I{;Uy61oix$ĿG Hk(& IVmrVL1ğ ;邦p cD{]\ZBJ9OlX>I7]8fb ozߐivX(߹ʄwͫKoy͊h]C m ~,E0*݌G}~ڛQza/QFiQ>ECb׮dĀVg__򹾼F<תT,: qFo߂$yiL.BV~+wMwYe4uxZ 7^ʼnʅ׎ĴUmd;778ރ)ROY›YdmA2j*jIֹ{nzNm.ڗl$Qhs\o %ӣ?@"1[!)3 H>re_t5xl9t⸁r@:M]uN6&8&BS*E}^r 9ZZGe/I9yd~mt "~5_UFSM( I€~<ⴅF 8>2D>Q7bc@^!mL!r񮳷z&9.oĻu؎&f?^C+tΊ6>ߝOZAڗ.w'u$cBYӓyo~/$^%(?(ku?x#Zd?ϛxاVH"ơ`H&,Gnѯ'؎:(,NFpn,q4`g2ĢF@;[?yi=pU]dFQ;kB,5%ڞc[RmHz:w&չ ~JMɹ9Jw >K3ĦDkA'CI /4w]| ᱐W6ho#h۽+ȋK{U jNp,΁6KT2&\C=N^7W8G gHmuxֶ`shN/?+CR;Y+Tr:>nR9)/@l*i@ Rktss_+FK] i&_ yH%X ~!; %Y7MU;7 niXdȢps6hTy=o\w ^Z9*CaBAaTZt\sBTtcdԩBED7dɕ۳Bfmv2-Ct ;GBV4_%% Y2ʏH0ևL{Ŝq6],;h|ɐw btHhNh5G; }Pt+ǣEw]K@6*#!+9Q/d@JA\*D]jA >L6 }Kh>d=\]S!UHC9xw8ђķGI|:|Y9 'U ˃^E Օf'{7#./BV yne`/_ ՜>|~`GbͰğ7|zZ>NJXo^فy{HfR@ d:2uA꤁JUp %9 s&BܫdrP-QE/BnIm'ٷT.d]뼈.sHS'd*/L.Vq˻:׈ʇyw-.?;EDCWeZZ{ hc3<> }H|Zx⟹we!EoZO*덫2Ur<*h&v?sChὠ: c[WoTrU"rODLEsodCvLD%ZzHt Y-/1FDEh_v)H+rB08d5%qtX'@66<ï3/}ܬ{{i]r{"({ɶ-ItB4 3#kv"GBmˈL<5&C,H//X>2 >7B}m=ZҗW>sӮ-_!>Hpgo?T,Qf2-/gFK_LFN0&_6 "Щtrŀ*|L_z@:oh}m7n9Asa([hQ_:!Z[ւK8#H8+C.Șf -Q`4CŧWyG*oGwRP4^7֓fuD6b5fUE"F#Fm3}J Sk1abP6=ԗOBwP6|EߗrF6nƑ(ة_@ p $)^|i4wwyý\Fcy(VQ(F){xh磬2: QP}[(UW59# #\>(z2$庫x(b<r;QBI;}V q椛i( 1Z|Fr ݠ^7Q8JPR)BAq^P*.OaX>G;U@ )ysK亸YͶ љ:,~ȧJ@#j%8ltXoU>hai8g<`dƭ6sk"G=k (mRށٍٶ(i盉C6Lg.Bم}Q@\?d6l9Co^$oM7QׯRLӬm9Vxy=U(^adg1E_{5E޵RL2نsZ^3Z~k$F~=t$YH9x}L+s^ODA]RY,|1RƣhSC.^|Z?&I*8< /^0~ de^znOAIBqкoG7FЋIr/OPV~k%TJ0袃>R8-.|tm`,h3|xMG@!FU_wcE5iچIivGZ%EQNe#'8<@ٴ2EW(~DX9t]ܘ!ylSakB r֑?JdlEEГlPOJ;j4E~q"v]{GݥQڧWrZuvY-"E+&h sm>Fk7P0ҩ(N7Q@ywxD GΑQ;Au7lsdІ}UOr-rB]y6g/Ikkһ? m:E&m*QB JP@(lk.םNţ ;A[pFyP,O(tzaȩ(z*P21Y/(2`3KUbA;B+/*(ZWclMgMP؇C忨/TnK} 3OBE~ߧB3-3Bn y;z?NGsƓ:^QE}n?xJff; b-w0tȍ~{ i zw Je?Dr֙{W?gPMF hq/@yOY40ض=LϵBүbse_1ky8 *%W$dmUSޠj;QaZr' dZfo^#uPuaʩ$(DeQLC -fPr-(Pd_;W{}zyP8ޱz+:HަyLEyèpd=v>jx'JiĢL CB/*vw{g5k*#JN@pq@uY0`Q ?{b8;En3ΔG^b Ww%@e}%TMV+;U?5OFIvkK=W'*w^Fů[m{f-QpgXn2*m6o]t6Oa[T}/Qׇ]J%xIUo̞fcBskX<;0m6݌SFO ѵ7лRxAZ*ϕ>+e+n;J;̐1jS#e0. Ioūy]yDT&5ωQBPyN5%/~=;Q6&8,*?Bԛe% ]30Azc!*N] _E먭2z~,Ng5: y<vF.Q*՝EC@և2ݼ;(.ʿ^ͷ1ʵ6N 7kv@}߽l0gy-F.JPߡҷ:cPO7: fp+sO'yEPcUHDGu`J-PGNfoqG'无.PZ-Cmx'>x4_fө YXкw++%5:Qq6;tN)ƥ3'CEOzUU.1c̓L~HA{ uǻs?*us"@] }_ >TR(vfwT0dCToTgs_{}vK﹦ $0RX`FP sݺ3P53lM"~įn8jAqR޷Qw_F.[l(-9; ? B%rII@9L{uIҊơv YqoԍzJ(U[' D1ל~XqHTK70XPL]J3mU\љŨ8|>,fS>B~TaE@;%DpN]vڰTW/dSBY~C>Tl4b\ϖW~@5'> 9ʾ{0t`O=l#*=mT8 LݽFEOφ~yoY:%)5?m5Kݞx< .(&p?|OoT 3M U)AJ?"~O L !2^F]jJ}֘StA`9+Pux?}'Q\&W䧛^H^CARP-1.0/DATA/.Rhistory0000755000076400007640000000057311263754702014375 0ustar volmelnvolmelnx = matrix(scan("x.dat"), ncol = 3, byrow = T) head(x) plot(x[,2], x[,3]) plot(x[,2], x[,3], cex = 0.5, pch = 19) x = matrix(scan("x.dat"), ncol = 3, byrow = T) plot(x[,2], x[,3], cex = 0.5, pch = 19) x = matrix(scan("x.dat"), ncol = 3, byrow = T) plot(x[,2], x[,3], cex = 0.5, pch = 19) x = matrix(scan("x.dat"), ncol = 3, byrow = T) plot(x[,2], x[,3], cex = 0.5, pch = 19) q() CARP-1.0/Matrix.c0000664000076400007640000001733611356774643013465 0ustar volmelnvolmeln /* CARP: Clustering Algorithms' Referee Package Copyright (C) 2010 Volodymyr Melnykov and Ranjan Maitra This program is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version. This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details. You should have received a copy of the GNU General Public License along with this program. If not, see . Authors' contact information: Volodymyr Melnykov Ranjan Maitra Volodymyr.Melnykov@ndsu.edu maitra@iastate.edu Department of Statistics Department of Statistics North Dakota State University Iowa State University Fargo, ND 58102 Ames, IA 50011 */ #include #include #include #include "array.h" #include "overlap.h" /* Multiplies matrices a and b and puts the result in c which should be pre-allocated */ void multiply(double **a, int arows, int acols, double **b, int brows, int bcols, double **c) { int i, j, k; for (i=0; i (*max)){ (*max) = x[i]; maxN = i; } } return maxN; } /* multiplies x by x' */ int vec11vecSQ(double *y, int p, double **Res){ int i,j; for (i=0;i. Authors' contact information: Volodymyr Melnykov Ranjan Maitra Volodymyr.Melnykov@ndsu.edu maitra@iastate.edu Department of Statistics Department of Statistics North Dakota State University Iowa State University Fargo, ND 58102 Ames, IA 50011 */ #include #include #include #include "array.h" #include "overlap.h" #define PI 3.141592653589793 #define EE 2.718281828459046 long int random(void); /* declare in order to avoid undef with ISO */ void srandom(unsigned int seed); /* declare in order to avoid undef with ISO */ double rxxx1; int ixxx1=0; /*************************************/ void setseed(unsigned int s) { srandom(s); } /*************************************/ long genseed(void) { return((long)random()); } /*************************************/ double runi(void) { return ((double)random() + 0.5)/((double)RAND_MAX + 1.0); } /*************************************/ double runir(double a,double b) { return (b-a)*runi()+a; } /*************************************/ double rnor(double mu,double sd) { double e,v1,v2,w; if(ixxx1==0){ do{ v1=2*runi()-1.; v2=2*runi()-1.; w=v1*v1+v2*v2; }while(w>1.); e=sqrt((-2.*log(w))/w); rxxx1=v1*e; ixxx1=1; return v2*e*sd+mu; } else{ ixxx1=0; return rxxx1*sd+mu; } } /*************************************/ double rexp(double lambda) { return -log(runi())/lambda; } /*************************************/ double rgamma(double alpha) { double r1,r2,aa,x,w,c1,c2,c3,c4,c5; if(alpha<=0.)return 0.; if(alpha == 1.)return rexp(1.); if(alpha<1){ aa=(alpha+EE)/EE; do{ r1=runi(); r2=runi(); if(r1>1./aa){ x = -log(aa*(1.-r1)/alpha); if(r22.5)r1=r2+c5*(1.-1.86*r1); }while(r1<=0 || r1 >= 1); w=c2*r2/r1; if(c3*r1+w+1/w <= c4)return c1*w; if(c3*log(r1)-log(w)+w<1)return c1*w; }while(r2<2); } exit(1); } /*************************************/ /* Generation of random realizations from Multinomial distribution * n - number of trials * K - dimensionality * pi - multinomial proportions * Nk - produced random realization */ void rmulti(int n, int K, double *pi, int *Nk){ int i, k; double X; double *cumPI; MAKE_VECTOR(cumPI, K); for(k=0; k < K; k++){ Nk[k] = 0; } cumPI[0] = pi[0]; for(k=1; k < K; k++){ cumPI[k] = cumPI[k-1] + pi[k]; } for (i = 0; i < n; i++){ X = runi(); for(k=0; k < K; k++){ if (X < cumPI[k]){ Nk[k] = Nk[k] + 1; break; } } } FREE_VECTOR(cumPI); } /*************************************/ /* Generation of random realizations from multivariate Gaussian distribution * n - sample size * K - dimensionality * MU - mean vector * S - covariance matrix * X - produced random realizations */ void rMVN(int n, int p, double *MU, double **S, double **X){ int i, j; double dtmt; double *Eig; double **Ga, **L, **Sh, **Z; MAKE_MATRIX(Sh, p, p); MAKE_VECTOR(Eig, p); MAKE_MATRIX(Ga, p, p); MAKE_MATRIX(L, p, p); MAKE_MATRIX(Z, n, p); cpy(S, p, p, Ga); cephes_symmeigens_down(p, Eig, Ga, &dtmt); Anull(L, p, p); for (i=0; i. Authors' contact information: Volodymyr Melnykov Ranjan Maitra Volodymyr.Melnykov@ndsu.edu maitra@iastate.edu Department of Statistics Department of Statistics North Dakota State University Iowa State University Fargo, ND 58102 Ames, IA 50011 */ #include #include #include #include #include #include #include "array.h" #include "overlap.h" int getopt(int argc, char *const argv[], const char *optstring); extern char *optarg; extern int optind, optopt; /* C-MixSim generates mixtures with prespecified level of overlap and datasets from these mixtures. */ int main(int argc, char *argv[]) { int i, k; int p, K, resN, method, lim, fail, sph, g, nf; double emax, PiLow, Ubound; double BarOmega, MaxOmega; int *rcMax; double *Pi, *pars; double **Mu, **OmegaMap; double ***S; int c, errflag; char **pp; pp = NULL; MAKE_VECTOR(pars, 2); /* DEFAULT PARAMETERS */ p = 2; K = 2; method = 0; /* 0 - average; 1 - maximum */ BarOmega = 0.0; MaxOmega = 0.0; sph = 0; /* 0 - nonspherical, 1 - spherical */ emax = 0.90; /* Max eccentricity */ PiLow = 1.0; /* Smallest mixing proportion */ Ubound = 1.0; /* Means are distributed by U(0,Ubound) */ resN = 100; /* Max # of dataset resimulations */ pars[0] = 0.0000001; /* eps */ pars[1] = 0.0000001; /* acc (Davies) */ lim = 1000000; /* lim (Davies) */ g = 0; /* sample size for generated datasets */ nf = 1; /* # of simulated mixtures */ char buff[200] = ""; char Path[100] = "DATA"; char PIfname[100] = "Pi.dat"; char MUfname[100] = "Mu.dat"; char Sfname[100] = "LTSigma.dat"; char overmap[100] = "overMap.dat"; char overbarmax[100] = "overBarMax.dat"; char dataX[100] = "x.dat"; char Nksizes[100] = "Nk.dat"; errflag = 0; /* READING PARAMETERS AND CHECKING THEM */ while ((c = getopt(argc, argv, ":b:m:p:K:se:z:u:r:a:l:P:M:S:D:I:X:W:C:n:#:")) != -1) { switch(c) { case 'b': /* average overlap */ BarOmega = strtod(optarg, pp); break; case 'm': /* maximum overlap */ MaxOmega = strtod(optarg, pp); break; case 'p': /* # of dimensions */ p = atoi(optarg); break; case 'K': /* # of components */ K = atoi(optarg); break; case 's': /* spherical cluster */ sph = 1; break; case 'e': /* maximum eccentricity */ emax = strtod(optarg, pp); break; case 'z': /* smallest mixing proportion */ PiLow = strtod(optarg, pp); break; case 'u': /* upper bound for Uniform */ Ubound = strtod(optarg, pp); break; case 'r': /* # of resimulations */ resN = atoi(optarg); break; case 'a': /*epsilon / accuracy */ pars[0] = strtod(optarg, pp); pars[1] = strtod(optarg, pp); break; case 'l': /* max # of integration terms */ lim = atoi(optarg); break; case 'P': /* Name of the file of mixing proportions */ strcpy(PIfname, optarg); break; case 'M': /* Name of the file with means */ strcpy(MUfname, optarg); break; case 'S': /* Name of the file with covariance matrices */ strcpy(Sfname, optarg); break; case 'D': /* Name of the working directory */ strcpy(Path, optarg); break; case 'I': /* Name of the file with true sample sizes */ strcpy(Nksizes, optarg); break; case 'X': /* Name of the output file */ strcpy(dataX, optarg); break; case 'W': /* Name of the file with overlap map */ strcpy(overmap, optarg); break; case 'C': /* Name of the file with overlap characteristics */ strcpy(overbarmax, optarg); break; case 'n': /* # of observations in simulated datasets */ g = atoi(optarg); break; case '#': /* # of simulated mixtures */ nf = atoi(optarg); break; case ':': /* option with no operand */ printf("Operand for option -%c is not specified correctly...\n", optopt); errflag = 1; break; case '?': printf("Unrecognized option: -%c\n", optopt); errflag = 2; } } /* DELAING WITH WRONG INPUT */ if (p < 1){ printf("Operand for option -p is not specified correctly...\n"); errflag = 1; } if (K < 1){ printf("Operand for option -K is not specified correctly...\n"); errflag = 1; } if ((sph != 0) & (sph != 1)){ printf("Operand for option -s is not specified correctly...\n"); errflag = 1; } if ((emax <= 0.0) | (emax >= 1.0)){ printf("Operand for option -e is not specified correctly...\n"); errflag = 1; } if ((PiLow <= 0.0) | (PiLow > 1.0)){ printf("Operand for option -z is not specified correctly...\n"); errflag = 1; } if (Ubound <= 0.0){ printf("Operand for option -u is not specified correctly...\n"); errflag = 1; } if (resN <= 0){ printf("Operand for option -r is not specified correctly...\n"); errflag = 1; } if (pars[0] <= 0.0){ printf("Operand for option -a is not specified correctly...\n"); errflag = 1; } if (lim <= 0){ printf("Operand for option -l is not specified correctly...\n"); errflag = 1; } if ((g < 0) & ((MaxOmega != 0.0) | (BarOmega != 0.0))){ printf("Operand for option -n is not specified correctly...\n"); errflag = 1; } if (nf < 1){ printf("Operand for option -# is not specified correctly...\n"); errflag = 1; } if (errflag != 0) exit(1); strcpy(buff, Path); strcat(buff, "/"); strcat(buff, PIfname); strcpy(PIfname, buff); strcpy(buff, Path); strcat(buff, "/"); strcat(buff, MUfname); strcpy(MUfname, buff); strcpy(buff, Path); strcat(buff, "/"); strcat(buff, Sfname); strcpy(Sfname, buff); strcpy(buff, Path); strcat(buff, "/"); strcat(buff, overmap); strcpy(overmap, buff); strcpy(buff, Path); strcat(buff, "/"); strcat(buff, overbarmax); strcpy(overbarmax, buff); strcpy(buff, Path); strcat(buff, "/"); strcat(buff, dataX); strcpy(dataX, buff); strcpy(buff, Path); strcat(buff, "/"); strcat(buff, Nksizes); strcpy(Nksizes, buff); srand(time(NULL)); MAKE_VECTOR(Pi, K); MAKE_MATRIX(Mu, K, p); MAKE_3ARRAY(S, K, p, p); MAKE_MATRIX(OmegaMap, K, K); MAKE_VECTOR(rcMax, 2); /* RUNNING A CORRESPONDING FUNCTION */ for (i=1; i <= nf; i++){ if (nf != 1){ printf("\n\nMIXTURE MODEL #%i:\n\n", i); } if (BarOmega == 0.0){ if (MaxOmega == 0.0){ /* AVERAGE AND MAXIMUM OVERLAPS ARE NOT SPECIFIED */ freadParameters(p, K, Pi, Mu, S, PIfname, MUfname, Sfname, i); ExactOverlap(p, K, Pi, Mu, S, pars, lim, OmegaMap, &BarOmega, &MaxOmega, rcMax); if ((BarOmega > MaxOmega) | (BarOmega*K*(K-1.0)/2 < MaxOmega) | (BarOmega < 0.0)){ printf("Incorrect parameters. "); printf("Check the options -p and -K...\n"); exit(2); } printOverlap(K, OmegaMap, BarOmega, MaxOmega, rcMax); fprintOverlap(K, OmegaMap, BarOmega, MaxOmega, rcMax, i, overmap, overbarmax); BarOmega = 0.0; MaxOmega = 0.0; if (g > 0){ double **Y; int *Nk; MAKE_MATRIX(Y, g, p); MAKE_VECTOR(Nk, K); /* dataset generation and printing parameters */ genData(p, K, Pi, Mu, S, g, Y, Nk); fprintData(g, p, K, Y, Nk, i, dataX, Nksizes); printf("Dataset with cluster sizes Nk = "); for (k=0; k. Authors' contact information: Volodymyr Melnykov Ranjan Maitra Volodymyr.Melnykov@ndsu.edu maitra@iastate.edu Department of Statistics Department of Statistics North Dakota State University Iowa State University Fargo, ND 58102 Ames, IA 50011 */ #include #include #include #include "array.h" /* fprintOverlap : writes overlaps to files * Parameters: * K - number of components * OmegaMap - map of misclassification probabilities * BarOmega - average overlap * MaxOmega - maximum overlap * rcMax - pair of components producing maximum overlap * fileN - number of the dataset * overmap - name of file with the map of misclassification probabilities * overbarmax - name of file with overlap characteristics */ void fprintOverlap(int K, double **OmegaMap, double BarOmega, double MaxOmega, int *rcMax, int fileN, char *overmap, char *overbarmax){ int i, j; FILE *fout; if (fileN == 1){ if (!(fout = fopen(overmap, "w"))){ fprintf(stderr, "Error: could not open file %s...\n", overmap); exit(EXIT_FAILURE); } } else { if (!(fout = fopen(overmap, "a"))){ fprintf(stderr, "Error: could not open file %s...\n", overmap); exit(EXIT_FAILURE); } } for (i=0; i tar -xzvf CARP_v1.0.tar.gz - compile files running the makefile: > make CARP To check the integrity of the package, run the command: > make check To remove the installed files, use the command: > make clean Usage: > ./CARP Parameters: -b: average overlap (no default value) -m: maximum overlap (no default value) -p: number of dimensions (2 by default) -K: number of mixing components (2 by default) -n: number of observations generated from every mixture (0 by default) -#: number of simulated mixtures (1 by default) -0: name of clustering program (if this option is not specified, only C-MixSim is run) -1: name of partition analyzing program (Adjusted Rand index is used by default; program AdjRand) -s: spherical covariance matrix structure (non-spherical by default if option is unspecified) -e: maximum eccentricity (0.90 by default) -z: smallest mixing proportion (equal mixing proportions 1/K by default) -u: upper bound for Uniform(0, upper-bound) distribution from which mean vectors are generated -r: maximum number of resimulations (100 by default) -a: accuracy of estimation (1e-06 by default) -l: maximum number of integration terms used by qfc function (1e06 by default) -P: file containing mixing proportions (Pi.dat by default) -M: file containing mean vectors (Mu.dat by default) -S: file containing covariance matrices in triangular form (LTSigma.dat by default) -D: working directory (DATA by default) -I: file containing numbers of observations generated from every cluster (Nk.dat by default) -i: file containing estimated classifications (idEst.dat by default) -X: file containing simulated datasets (x.dat by default) -W: file containing maps of pairwise overlaps (overMap.dat by default) -C: file containing characteristics of simulated mixtures (overBarMax.dat by default) -R: file containing index values (Adjusted Rand index and AR.dat by default) Details: Upon launching CARP, three stages of the program are processed. The first stage is reponsible for the simulation of Gaussian mixtures with prespecified level of complexity expressed in terms of pairwise overlap and for the generation of datasets from these mixtures. At this stage, CARP invokes C-MixSim. The admissible options are the following: -b, -m, -p, -K, -n, -#, -s, -e, -z, -u, -r, -a, -l, -P, -M, -S, -D, -I, -X, -W, -C. At the second stage, a user-specified clustering method has to be run for the simulated datasets. For illustration purposes, the hierarchical clustering algorithm hierclust is employed. User’s clustering priogram has to be able to accept the following options: -p, -K, -n, -#, -D, -i, -X. The program has to obtain partitionings and write them into the file specified by the option -i. At the third stage, true classification vectors are compared with obtained ones to assess the performance of the clustering algorithm. The Adjusted Rand index is incorporated as a default measure of similarity between the two classifications. A user, however, can specify some other program for comparing the obtained and true partitionings. Then, the program should comply with the following options: -K, -n, -#, -D, -I, -i, -R. Calculated values of the provided similarity measure have to be written to the file specified by the option -R. If both options -b and -m are specified, CARP produces a mixture satisfy ing both characteristics, average and maximum overlap. If one option, -b or -m is specified, a mixture satisfying the prespecified value is generated. The working directory is specified by the option -D. The obtained parameters will be saved to the files specified by options -P, -M, and -S while samples drawn from the mixtures will be saved into the file specified by the option -X. The file provided with the option -I contains sample sizes for every mixture component. In addition, the map of misclassification probabilities will be saved into the file specified by the option -W while average and maximum overlaps as well as the row and column numbers of the components that produced maximum overlap are stored in the file given by the option -C. The element with the index (i, j) in the misclassification map represents the probability that X simulated from the ith component is classified to the jth component. If both options, -b and -m are not provided, CARP reads parameters from the files specified by the options -P, -M, -S and computes misclassification probabilities for the components of the given mixture model. Note that the options -p and -K have to be correctly specified. The options -n and -# specify the sample size of a sample drawn from every generated mixture and the number of such mixtures correspondingly. If it is desired to use datasets stored in a file specified by the option -X, the sample size should be given as a negative number, for example: -n-100 instead of the traditional -n100. In that case, new datasets will not be simulated. Note that this capability is available only in the mode when mixture parameters are not simulated but rather read from files. This capability can be useful in the case when it is desired to run several clustering algorithms on the same set of simulated data. Another application of CARP’s capabilities is for assessing the level of clustering complexity of derived and existing groups in datasets based on the set of estimated mixture model parameters. CARP also has a capability to simulate finite mixture models and datasets from them without launching the second and third stages. Thus, the simulation stage would be the only one to run. In order to do so, the option -0 has to be omitted. Simple Examples: % simulate a 3-dimensional 4-component mixture with spherical covariance matrices, equal mixing proportions, maximum overlap 0.1 and average overlap 0.05; generate a sample of size 100 and analyze using hierclust > ./CARP -m0.1 -b0.05 -p3 -K4 -s -n100 -0hierclust % simulate a 2-dimensional 4-component mixture with maximum overlap 0.1 and mixing proportions greater or equal than 0.10; generate a sample of size 100 > ./CARP -m0.1 -p2 -K4 -z0.10 -n100 % compute the overlap for the 2-dimensional 4-component mixture with parameters specified in DATA/Pi.dat, DATA/Mu.dat and DATA/LTSigma.dat; generate a sample of size 100 > ./CARP -p2 -K4 -n100 % get help > ./CARP See the manual for more comprehensive examples. References: L. Hubert and P. Arabie. Comparing partitions. Journal of Classification, 2:193–218, 1985. R. Maitra and V. Melnykov. Simulating data to study performance of finite mixture modeling and clustering algorithms. Accepted for publication in the Journal of Computational and Graphical Statistics, 2010. doi: 10.1198/jcgs.2009.08054.